Breeding Bird Surveys - Cleveland Metroparks
Nathan W. Byer
5/16/2022
The dataset
Tim Krynak, Erik Shaffer, Jennifer Brumsfield, and other Cleveland Metroparks and affiliated colleagues conducted Breeding Bird Surveys (BBS) at a total of 211 routes between 2017 and 2021. 164 of these routes overlapped with long-term vegetation monitoring efforts through the Cleveland Metroparks Plant Community Assessment Program (PCAP). This provided an excellent opportunity to explore the potential relationships between vegetative composition, landscape structure, and breeding bird ecology and community composition.
Initial processing
First, we just read in the full BBS dataset. For our purposes, we will remove any BBS observations not paired with PCAP plots. We can also visualize the distribution of survey effort through time.
CMP_bbs_pcap_join<-read.csv("CMP_bbs_pcapjoin.csv",header=T)
CMP_bbs_pcap_join_pcaponly<-CMP_bbs_pcap_join[!is.na(CMP_bbs_pcap_join$Plot.ID),]
datetime<-strptime(CMP_bbs_pcap_join_pcaponly$survey_cre,format="%m/%d/%Y %H:%M")
hist(datetime,breaks="months")
Next, we need to start grouping this dataset for downsteram analysis. This involves summarizing counts of each species for each plot. We also replace any NA values in the data frame with zeroes - this will be important for downstream applications, as we want to distinguish between unobserved species at a plot and check (typically coded 0s) and checks that were “skipped” (NAs).
CMP_bbs_pcap_join_pcaponly_pcapsum<-CMP_bbs_pcap_join_pcaponly %>%
group_by(simp_plot, reservatio,species_4_) %>%
summarise(n=n(),longitude=mean(xcoord), latitude=mean(ycoord),vibifq=mean(vibi_fq),vibifqnotrees=mean(vibi_fq_no_trees),carex=mean(carex_metric_value),cyper=mean(cyperaceae_metric_value),
dicot=mean(dicot_metric_value),shade=mean(shade_metric_value),shrub=mean(shrub_metric_value),
hydro=mean(hydrophyte_metric_value),svp=mean(svp_metric_value),ap=mean(ap_ratio_metric_value),
fqi=mean(fqai_score),bryo=mean(bryophyte_metric_value),hydro_rel=mean(hydrophyte_rel_cov_metric_value),
hydro_rel_notrees=mean(hydrophyte_rel_cov_no_trees),
sens_rel=mean(sensitive_rel_cov_metric_value),
sens_rel_notree=mean(sensitive_rel_cov_no_trees),
tol_rel=mean(tolerant_rel_cov_metric_value),
tol_rel_notree=mean(tolerant_rel_cov_no_trees),
inv=mean(invasive_graminoids_metric_value),
small_tree=mean(small_tree_metric_value),
subcanopy=mean(subcanopy_iv),
canopy=mean(canopy_iv.1)
) %>%
spread(species_4_,n) %>%
replace(is.na(.),0)
CMP_bbs_pcap_join_pcaponly_pcapsum<-column_to_rownames(CMP_bbs_pcap_join_pcaponly_pcapsum,var="simp_plot")We can also take a look at this dataset to ensure that it is what we expect.
head(CMP_bbs_pcap_join_pcaponly_pcapsum)## reservatio longitude latitude vibifq vibifqnotrees carex cyper dicot shade shrub hydro svp ap fqi bryo hydro_rel hydro_rel_notrees sens_rel sens_rel_notree tol_rel tol_rel_notree inv small_tree subcanopy canopy V1 ACFL ALFL AMCR AMGO AMKE AMRE AMRO AMWO BADO
## 1002 HI 2181968 561930.7 47.26427 24.35064 9 9 30 25 1 9 1 0.05882353 19.47063 0.005272844 0.01581853 0.03168627 0.03962464 0.000313725 0.41913886 0.76606536 0 0.13936430 0.03419060 0.1922540 0 1 0 1 2 0 2 0 0 0
## 1003 MS 2150599 613095.3 76.16889 52.67227 4 4 22 22 1 2 1 0.13333333 22.98097 0.009033832 0.01806766 0.06677053 0.48845928 0.051747162 0.11096557 0.04284442 0 0.04845815 0.03497367 0.2054565 1 5 0 0 4 0 0 4 0 0
## 1004 BE 2227458 626184.1 37.88422 35.10323 2 2 21 23 1 9 4 0.19047619 15.08659 0.039605806 0.42634485 0.55871039 0.01397852 0.009159187 0.19872796 0.26042621 0 0.10869565 0.46386133 0.2174377 1 5 0 3 1 0 0 5 0 0
## 1005 HI 2191849 580242.5 82.91987 65.84743 11 11 47 52 2 18 3 0.08000000 30.51229 0.007461016 0.06013579 0.12854181 0.16854436 0.088512041 0.13320401 0.19685290 0 0.11436170 0.10777783 0.1757489 0 3 0 2 2 0 0 2 0 0
## 1006 BW 2116879 635151.5 72.07176 64.74055 3 3 21 23 1 6 2 0.06666667 26.58228 0.008306874 0.11007532 0.25470932 0.25985749 0.046004015 0.12386473 0.06877109 0 0.13235294 0.22866471 0.1707812 0 2 0 1 1 0 0 4 0 0
## 1008 WC 2190113 621267.9 56.30316 25.13574 4 4 11 13 0 3 0 0.25000000 16.84286 0.010954302 0.05477151 0.25875453 0.32594526 0.051750906 0.06057729 0.28618251 0 0.13084112 0.04519758 0.2905217 0 0 0 1 3 0 0 3 0 0
## BAEA BANS BAOR BARS BCCH BEKI BGGN BHCO BHVI BLJA BOBO BRCR BRTH BTNW BWHA BWWA CANG CARW CEDW CERW CHSP CHSW CLSW COGR COHA COYE DEJU DOWO EABL EAKI EAME EAPH EATO EAWP EUST FISP GBHE GCFL GHOW GRCA GRHE HAWO HERG HETH HOFI HOSP HOWA HOWR INBU KILL LOWA MALL MODO NOCA NOFL NOPA NRWS OROR OVEN
## 1002 0 0 0 0 1 0 1 2 0 1 0 0 0 0 0 1 0 0 1 0 0 1 0 0 0 0 0 2 0 0 0 0 3 2 0 0 0 0 0 0 0 0 0 0 0 0 3 0 2 0 0 0 0 3 0 0 0 0 0
## 1003 0 0 2 0 2 0 1 2 0 2 0 0 0 0 0 0 0 0 1 0 0 0 0 2 0 0 0 1 0 0 0 0 0 3 0 0 0 1 0 0 0 0 0 0 0 0 0 1 0 0 0 0 0 5 0 0 0 0 0
## 1004 0 0 0 0 6 1 2 3 0 6 0 0 0 0 0 0 1 0 0 0 0 2 0 0 0 0 1 3 0 0 0 0 0 7 0 0 0 1 0 0 0 1 0 0 0 0 0 1 2 0 0 0 0 3 1 0 0 0 0
## 1005 0 0 1 0 1 0 0 2 0 3 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 2 0 2 0 1 0 0 1 3 0 0 0 1 0 3 0 0 0 0 0 0 1 0 0 0 0 0 0 3 0 0 0 0 0
## 1006 0 0 1 0 6 0 2 0 0 7 0 0 0 0 0 0 0 0 0 0 0 1 0 1 0 0 0 5 0 0 0 0 0 7 0 0 0 3 0 0 0 3 0 0 0 0 0 0 0 0 0 0 0 5 1 0 0 0 0
## 1008 0 0 0 0 1 0 1 0 0 2 0 0 0 0 0 0 0 0 1 0 0 0 0 2 0 0 0 3 0 0 0 0 0 3 1 0 0 1 0 1 0 0 0 0 1 0 0 3 0 0 0 0 0 2 1 0 0 0 0
## PIWO PUFI PUMA RBGR RBGU RBNU RBWO REVI RHWO ROPI RSHA RTHA RTHU RWBL SAVS SCTA SOSP SWSP TEWA TRES TUTI TUVU VEER VIRA WAVI WBNU WEVI WIFL WITU WIWR WODU WOTH YBCU YEWA YTVI YTWA
## 1002 0 0 0 0 0 0 0 2 0 0 0 0 0 0 0 0 0 0 0 0 1 0 0 0 0 1 0 0 0 0 0 2 0 0 1 0
## 1003 0 0 0 0 0 0 1 5 0 0 0 0 0 0 0 2 0 0 0 0 2 0 1 0 0 5 0 0 0 0 0 4 0 0 0 0
## 1004 0 0 0 0 2 0 5 7 0 0 1 0 0 0 0 5 0 0 0 0 6 0 0 0 0 4 0 0 0 0 0 3 0 0 0 0
## 1005 1 0 0 1 0 0 3 1 0 0 0 0 1 3 0 1 2 3 0 0 3 0 0 0 0 3 0 1 0 0 0 2 0 0 0 0
## 1006 4 0 0 4 0 0 6 6 0 0 0 0 0 0 0 6 0 0 0 0 4 0 2 0 0 3 0 0 0 1 0 6 0 0 0 0
## 1008 0 0 0 0 0 0 2 2 2 0 0 1 0 0 0 0 0 0 0 0 3 0 0 0 0 2 0 0 0 0 0 0 0 0 2 0This is a pretty unwieldly dataset at present - we have 164 rows (representing our BBS/PCAP plots) and 131 columns. Out of those 131 columns, 104 are unique species (designated with four-letter codes), so most of the first 27 covariates are essentially our “predictors” - in this case, our PCAP data.
Let’s do a bit more processing on this dataset to get it in shape for downstream analyses. First, we will separate out our bird community data.
birdcomm<-CMP_bbs_pcap_join_pcaponly_pcapsum[,27:(ncol(CMP_bbs_pcap_join_pcaponly_pcapsum)-1)]
birdcomm[is.na(birdcomm)] <- 0Next, we separate out our PCAP data.
vegecomm<-CMP_bbs_pcap_join_pcaponly_pcapsum[,4:25]Since we are using the PCAP dataset as predictors here, we should likely check for correlations between these covariates. I have selected a (semi-)arbitrary cut-off of 0.6 to identify covariates as correlated; we then drop any covariates that are likely redundant.
vegecomm_cor<-cor(vegecomm)
set.seed(12345)
corvar<-findCorrelation(vegecomm_cor,cutoff=0.6)
vegecomm_drop<-vegecomm[,-corvar]Now, we have:
- a vegecomm dataset - which representing PCAP measurements for each of the 164 BBS/PCAP plots.
- a birdcomm dataset - representing presence/absence of each bird species for each survey in each plot.
Note that we are not explicitly considering differences in bird counts (or occupancy, or detectability) through time at this point - we will get to that later. For now, we just want to explore broad associations between these two datasets!
Note as of 6/13/2022 - we will also want to only work with only binary data.
Community analyses
When working with datasets like these, we first have to acknowledge that these datasets are multivariate - with many possible predictor covariates (for example, plot-level characteristics) and many possible response covariates (for example, presence/absence of species in a plot). Fortunately, community ecology as a field has developed a number of standard approaches for exploring multivariate datasets, of which the most common is (arguably!) ordination. In essence, ordination techniques try to shrink down highly multidimensional datasets into a more manageable number of dimensions - often just a handful (and for initial plotting, typically just two). There are many, many possible ordination tools available for exploring community ecological datasets, and we will explore a few.
unconstrained ordination
In this section, we will be exploring ways of exploring the bird dataset on its own - basically, without any predictors.
First, we must reclassify our bird community data onto a binary scale - ranging from 0 (no detection) to 1 (detection).
birdcomm_bin<-birdcomm
birdcomm_bin<- birdcomm_bin %>% mutate_if(is.numeric, ~1 * (. >= 1))One of the most common approaches for unconstrained ordination is called Principal Components Analysis (PCA). I have illustrated two approaches for this below - using either the prcomp() function or vegan’s rda() function (which we will use a bit more below!).
PCA.res<-prcomp(birdcomm_bin, scale=T)
biplot(PCA.res)
title("PCA (prcomp version)")
Personally, I find rda() to be much easier to use, so I tend to use it preferentially. Here is a slightly more complicated plot compared to the above - this time with convex hulls and points colored by reservation.
reserv<-factor(toupper(CMP_bbs_pcap_join_pcaponly_pcapsum$reservatio))
colpal_reserv<-rainbow(length(levels(reserv)))
p_unconc<-rda(birdcomm_bin,scale=T)
PCAscores <- scores(p_unconc, display = "sites") %>%
as.data.frame() %>%
rownames_to_column("site")
PCAvect <- scores(p_unconc, display = "species") %>%
as.data.frame()
reservplot<-ordiplot(p_unconc,scaling=3,type="n")
points(reservplot, "sites",pch=21, col=colpal_reserv[reserv])
ordihull(reservplot,groups = reserv, col=colpal_reserv, lwd=3)
title("PCA (vegan rda version)")
legend(x="right", legend=levels(reserv), col=colpal_reserv,pch = 1,cex=0.5) 
while a bit hard to interpret, plots within Rocky River (RR) seem to be the most variable - and it might be that we see higher diversity in that area. We will see as we proceed if that is true.
metric multidimensional scaling (Principal Coordinates Analysis)
You can think of this as a distance-based analog to PCA. In this case, since we are using Bray-Curtis distances on binary data (presence/absence), this is equivalent to Sorensen’s index.
non-metric multidimensional scaling
bird_nmds<-metaMDS(birdcomm_bin,distance = "bray",trymax = 400,k=3)## Run 0 stress 0.1751474
## Run 1 stress 0.1761345
## Run 2 stress 0.1732776
## ... New best solution
## ... Procrustes: rmse 0.03981098 max resid 0.1508202
## Run 3 stress 0.1743699
## Run 4 stress 0.1738508
## Run 5 stress 0.1725355
## ... New best solution
## ... Procrustes: rmse 0.01742725 max resid 0.1440867
## Run 6 stress 0.1733826
## Run 7 stress 0.1753801
## Run 8 stress 0.1727212
## ... Procrustes: rmse 0.02129917 max resid 0.1476814
## Run 9 stress 0.172957
## ... Procrustes: rmse 0.02323107 max resid 0.1467514
## Run 10 stress 0.1724776
## ... New best solution
## ... Procrustes: rmse 0.00256073 max resid 0.01515115
## Run 11 stress 0.1725733
## ... Procrustes: rmse 0.00451041 max resid 0.03826016
## Run 12 stress 0.172205
## ... New best solution
## ... Procrustes: rmse 0.01955409 max resid 0.1465464
## Run 13 stress 0.1732051
## Run 14 stress 0.175053
## Run 15 stress 0.1752042
## Run 16 stress 0.1763912
## Run 17 stress 0.1745713
## Run 18 stress 0.1729667
## Run 19 stress 0.1734391
## Run 20 stress 0.1729283
## Run 21 stress 0.1732458
## Run 22 stress 0.1752081
## Run 23 stress 0.1754885
## Run 24 stress 0.1738273
## Run 25 stress 0.1721794
## ... New best solution
## ... Procrustes: rmse 0.003391018 max resid 0.03187639
## Run 26 stress 0.173727
## Run 27 stress 0.1755421
## Run 28 stress 0.172803
## Run 29 stress 0.1728229
## Run 30 stress 0.1725033
## ... Procrustes: rmse 0.02094271 max resid 0.1538802
## Run 31 stress 0.173439
## Run 32 stress 0.1751567
## Run 33 stress 0.1724449
## ... Procrustes: rmse 0.01094501 max resid 0.05449719
## Run 34 stress 0.1722481
## ... Procrustes: rmse 0.0153566 max resid 0.1417752
## Run 35 stress 0.1727234
## Run 36 stress 0.1748387
## Run 37 stress 0.1735728
## Run 38 stress 0.1752443
## Run 39 stress 0.1733246
## Run 40 stress 0.1726289
## ... Procrustes: rmse 0.02159124 max resid 0.1466149
## Run 41 stress 0.1729111
## Run 42 stress 0.173677
## Run 43 stress 0.1747492
## Run 44 stress 0.1721845
## ... Procrustes: rmse 0.004780107 max resid 0.04519247
## Run 45 stress 0.1744348
## Run 46 stress 0.1728829
## Run 47 stress 0.1745146
## Run 48 stress 0.1725974
## ... Procrustes: rmse 0.008605116 max resid 0.05020094
## Run 49 stress 0.1767774
## Run 50 stress 0.1738603
## Run 51 stress 0.1725219
## ... Procrustes: rmse 0.01557841 max resid 0.1479282
## Run 52 stress 0.173304
## Run 53 stress 0.1724993
## ... Procrustes: rmse 0.01543618 max resid 0.145884
## Run 54 stress 0.1764419
## Run 55 stress 0.1734277
## Run 56 stress 0.1731943
## Run 57 stress 0.1725213
## ... Procrustes: rmse 0.01543772 max resid 0.1475148
## Run 58 stress 0.1784585
## Run 59 stress 0.1740317
## Run 60 stress 0.172172
## ... New best solution
## ... Procrustes: rmse 0.002653639 max resid 0.02017705
## Run 61 stress 0.1759132
## Run 62 stress 0.1729089
## Run 63 stress 0.173989
## Run 64 stress 0.1758112
## Run 65 stress 0.1736825
## Run 66 stress 0.1733895
## Run 67 stress 0.17459
## Run 68 stress 0.1728892
## Run 69 stress 0.1736083
## Run 70 stress 0.1750135
## Run 71 stress 0.1750304
## Run 72 stress 0.1725164
## ... Procrustes: rmse 0.01549591 max resid 0.1468597
## Run 73 stress 0.1737165
## Run 74 stress 0.1742695
## Run 75 stress 0.1736979
## Run 76 stress 0.1724548
## ... Procrustes: rmse 0.01669034 max resid 0.1456585
## Run 77 stress 0.1728158
## Run 78 stress 0.1743856
## Run 79 stress 0.1754929
## Run 80 stress 0.1726014
## ... Procrustes: rmse 0.01458213 max resid 0.1469235
## Run 81 stress 0.1729228
## Run 82 stress 0.1724446
## ... Procrustes: rmse 0.01699636 max resid 0.1432058
## Run 83 stress 0.1738788
## Run 84 stress 0.1733452
## Run 85 stress 0.1733887
## Run 86 stress 0.1741618
## Run 87 stress 0.1749985
## Run 88 stress 0.1753223
## Run 89 stress 0.1772594
## Run 90 stress 0.1743855
## Run 91 stress 0.1755166
## Run 92 stress 0.1759072
## Run 93 stress 0.1752903
## Run 94 stress 0.1725475
## ... Procrustes: rmse 0.01843695 max resid 0.1509723
## Run 95 stress 0.1730089
## Run 96 stress 0.1743684
## Run 97 stress 0.174357
## Run 98 stress 0.1722625
## ... Procrustes: rmse 0.01545736 max resid 0.1381
## Run 99 stress 0.172896
## Run 100 stress 0.1732551
## Run 101 stress 0.1757211
## Run 102 stress 0.1750701
## Run 103 stress 0.1728711
## Run 104 stress 0.1759084
## Run 105 stress 0.1727175
## Run 106 stress 0.1756873
## Run 107 stress 0.1738645
## Run 108 stress 0.1742867
## Run 109 stress 0.1740164
## Run 110 stress 0.1752175
## Run 111 stress 0.1723468
## ... Procrustes: rmse 0.01984014 max resid 0.1459213
## Run 112 stress 0.1760383
## Run 113 stress 0.1728909
## Run 114 stress 0.1725117
## ... Procrustes: rmse 0.01563108 max resid 0.1470691
## Run 115 stress 0.1744335
## Run 116 stress 0.1752269
## Run 117 stress 0.1736964
## Run 118 stress 0.175521
## Run 119 stress 0.1727817
## Run 120 stress 0.1727144
## Run 121 stress 0.172672
## ... Procrustes: rmse 0.02113729 max resid 0.1400382
## Run 122 stress 0.1765034
## Run 123 stress 0.1752396
## Run 124 stress 0.1727298
## Run 125 stress 0.1740477
## Run 126 stress 0.172341
## ... Procrustes: rmse 0.01964859 max resid 0.145326
## Run 127 stress 0.174959
## Run 128 stress 0.1732102
## Run 129 stress 0.172554
## ... Procrustes: rmse 0.01377573 max resid 0.1275104
## Run 130 stress 0.1725232
## ... Procrustes: rmse 0.01199473 max resid 0.1173779
## Run 131 stress 0.1757938
## Run 132 stress 0.1769059
## Run 133 stress 0.1725009
## ... Procrustes: rmse 0.01591585 max resid 0.1463755
## Run 134 stress 0.1745496
## Run 135 stress 0.172888
## Run 136 stress 0.1762493
## Run 137 stress 0.1735169
## Run 138 stress 0.1725253
## ... Procrustes: rmse 0.01301501 max resid 0.1177242
## Run 139 stress 0.1755574
## Run 140 stress 0.1748397
## Run 141 stress 0.1741356
## Run 142 stress 0.1765357
## Run 143 stress 0.1731915
## Run 144 stress 0.172731
## Run 145 stress 0.1721813
## ... Procrustes: rmse 0.002955712 max resid 0.02239333
## Run 146 stress 0.1723639
## ... Procrustes: rmse 0.02035249 max resid 0.1485191
## Run 147 stress 0.1721948
## ... Procrustes: rmse 0.005701654 max resid 0.04621034
## Run 148 stress 0.1734274
## Run 149 stress 0.1728666
## Run 150 stress 0.1730735
## Run 151 stress 0.1739775
## Run 152 stress 0.1731944
## Run 153 stress 0.1722683
## ... Procrustes: rmse 0.01543467 max resid 0.1331563
## Run 154 stress 0.173358
## Run 155 stress 0.1736016
## Run 156 stress 0.1761904
## Run 157 stress 0.1725567
## ... Procrustes: rmse 0.01479097 max resid 0.1455197
## Run 158 stress 0.1746467
## Run 159 stress 0.1744006
## Run 160 stress 0.1742816
## Run 161 stress 0.1724327
## ... Procrustes: rmse 0.01706571 max resid 0.1450522
## Run 162 stress 0.1748837
## Run 163 stress 0.1743231
## Run 164 stress 0.1736804
## Run 165 stress 0.1753484
## Run 166 stress 0.1738681
## Run 167 stress 0.1735758
## Run 168 stress 0.1741785
## Run 169 stress 0.1755309
## Run 170 stress 0.1735338
## Run 171 stress 0.1722336
## ... Procrustes: rmse 0.01601707 max resid 0.1399591
## Run 172 stress 0.1723744
## ... Procrustes: rmse 0.02093664 max resid 0.1545463
## Run 173 stress 0.1723332
## ... Procrustes: rmse 0.0195882 max resid 0.14724
## Run 174 stress 0.175532
## Run 175 stress 0.1727458
## Run 176 stress 0.1732635
## Run 177 stress 0.1759102
## Run 178 stress 0.1740541
## Run 179 stress 0.1721916
## ... Procrustes: rmse 0.005767637 max resid 0.04636342
## Run 180 stress 0.1749051
## Run 181 stress 0.173469
## Run 182 stress 0.1726598
## ... Procrustes: rmse 0.02291106 max resid 0.1514561
## Run 183 stress 0.1724743
## ... Procrustes: rmse 0.0203851 max resid 0.1472383
## Run 184 stress 0.1752857
## Run 185 stress 0.1739626
## Run 186 stress 0.1738131
## Run 187 stress 0.1721772
## ... Procrustes: rmse 0.001878673 max resid 0.01768344
## Run 188 stress 0.1745026
## Run 189 stress 0.1726038
## ... Procrustes: rmse 0.02171574 max resid 0.1486159
## Run 190 stress 0.1742865
## Run 191 stress 0.1722344
## ... Procrustes: rmse 0.01624271 max resid 0.1427985
## Run 192 stress 0.1743934
## Run 193 stress 0.1725121
## ... Procrustes: rmse 0.01557617 max resid 0.1470384
## Run 194 stress 0.1740988
## Run 195 stress 0.1729496
## Run 196 stress 0.172682
## Run 197 stress 0.1746104
## Run 198 stress 0.1758245
## Run 199 stress 0.1762314
## Run 200 stress 0.1729143
## Run 201 stress 0.1745287
## Run 202 stress 0.1725161
## ... Procrustes: rmse 0.01567508 max resid 0.1477249
## Run 203 stress 0.1751867
## Run 204 stress 0.1733611
## Run 205 stress 0.172937
## Run 206 stress 0.1740611
## Run 207 stress 0.1722337
## ... Procrustes: rmse 0.01782391 max resid 0.1466053
## Run 208 stress 0.1745937
## Run 209 stress 0.173274
## Run 210 stress 0.1748937
## Run 211 stress 0.1739471
## Run 212 stress 0.1734346
## Run 213 stress 0.1742989
## Run 214 stress 0.1733951
## Run 215 stress 0.1731381
## Run 216 stress 0.1722191
## ... Procrustes: rmse 0.006453694 max resid 0.04636308
## Run 217 stress 0.1744902
## Run 218 stress 0.1755554
## Run 219 stress 0.1726841
## Run 220 stress 0.1727686
## Run 221 stress 0.1722391
## ... Procrustes: rmse 0.01774966 max resid 0.147681
## Run 222 stress 0.1767903
## Run 223 stress 0.1752083
## Run 224 stress 0.1726425
## ... Procrustes: rmse 0.01806157 max resid 0.1415587
## Run 225 stress 0.174496
## Run 226 stress 0.1729073
## Run 227 stress 0.1756222
## Run 228 stress 0.174881
## Run 229 stress 0.1736295
## Run 230 stress 0.1749679
## Run 231 stress 0.1725654
## ... Procrustes: rmse 0.01652435 max resid 0.1563714
## Run 232 stress 0.1749454
## Run 233 stress 0.175765
## Run 234 stress 0.1723714
## ... Procrustes: rmse 0.01981196 max resid 0.144888
## Run 235 stress 0.1751377
## Run 236 stress 0.1739592
## Run 237 stress 0.1723483
## ... Procrustes: rmse 0.01780819 max resid 0.1386382
## Run 238 stress 0.1725032
## ... Procrustes: rmse 0.01570548 max resid 0.1460301
## Run 239 stress 0.1739827
## Run 240 stress 0.1763483
## Run 241 stress 0.1730713
## Run 242 stress 0.1733056
## Run 243 stress 0.1727
## Run 244 stress 0.175973
## Run 245 stress 0.1725987
## ... Procrustes: rmse 0.01640206 max resid 0.1588538
## Run 246 stress 0.1730129
## Run 247 stress 0.175144
## Run 248 stress 0.1779461
## Run 249 stress 0.1783918
## Run 250 stress 0.1729347
## Run 251 stress 0.1743098
## Run 252 stress 0.1723768
## ... Procrustes: rmse 0.02059765 max resid 0.1504712
## Run 253 stress 0.1733464
## Run 254 stress 0.1748556
## Run 255 stress 0.1733471
## Run 256 stress 0.1754993
## Run 257 stress 0.1759528
## Run 258 stress 0.1753172
## Run 259 stress 0.1747468
## Run 260 stress 0.176772
## Run 261 stress 0.1746528
## Run 262 stress 0.173722
## Run 263 stress 0.1721774
## ... Procrustes: rmse 0.002415078 max resid 0.01975738
## Run 264 stress 0.1724053
## ... Procrustes: rmse 0.02043148 max resid 0.1515156
## Run 265 stress 0.1738366
## Run 266 stress 0.1754295
## Run 267 stress 0.1729009
## Run 268 stress 0.1729242
## Run 269 stress 0.1775915
## Run 270 stress 0.1727292
## Run 271 stress 0.1728591
## Run 272 stress 0.1759786
## Run 273 stress 0.1724741
## ... Procrustes: rmse 0.02040129 max resid 0.1473668
## Run 274 stress 0.1733485
## Run 275 stress 0.1736517
## Run 276 stress 0.1741278
## Run 277 stress 0.1721775
## ... Procrustes: rmse 0.002352542 max resid 0.02004201
## Run 278 stress 0.1742175
## Run 279 stress 0.1725201
## ... Procrustes: rmse 0.01383473 max resid 0.1180847
## Run 280 stress 0.1727011
## Run 281 stress 0.1725257
## ... Procrustes: rmse 0.020634 max resid 0.1548589
## Run 282 stress 0.174761
## Run 283 stress 0.1754632
## Run 284 stress 0.1751458
## Run 285 stress 0.1724602
## ... Procrustes: rmse 0.02056143 max resid 0.1475234
## Run 286 stress 0.1734248
## Run 287 stress 0.1726185
## ... Procrustes: rmse 0.02165034 max resid 0.141876
## Run 288 stress 0.1728716
## Run 289 stress 0.1745318
## Run 290 stress 0.1737329
## Run 291 stress 0.1753969
## Run 292 stress 0.1760232
## Run 293 stress 0.1723744
## ... Procrustes: rmse 0.0056035 max resid 0.03506882
## Run 294 stress 0.1731454
## Run 295 stress 0.1721988
## ... Procrustes: rmse 0.004916334 max resid 0.04590163
## Run 296 stress 0.173888
## Run 297 stress 0.1763478
## Run 298 stress 0.1722104
## ... Procrustes: rmse 0.004017098 max resid 0.0320913
## Run 299 stress 0.1725101
## ... Procrustes: rmse 0.008473902 max resid 0.05683405
## Run 300 stress 0.1755023
## Run 301 stress 0.1726061
## ... Procrustes: rmse 0.01799632 max resid 0.1455019
## Run 302 stress 0.1727611
## Run 303 stress 0.1736767
## Run 304 stress 0.1728906
## Run 305 stress 0.1740513
## Run 306 stress 0.1735696
## Run 307 stress 0.1739728
## Run 308 stress 0.1787724
## Run 309 stress 0.1736607
## Run 310 stress 0.1737475
## Run 311 stress 0.1727552
## Run 312 stress 0.1726611
## ... Procrustes: rmse 0.02117768 max resid 0.1404394
## Run 313 stress 0.1733226
## Run 314 stress 0.1727581
## Run 315 stress 0.1732464
## Run 316 stress 0.1723731
## ... Procrustes: rmse 0.02080382 max resid 0.151506
## Run 317 stress 0.1752211
## Run 318 stress 0.175093
## Run 319 stress 0.1766787
## Run 320 stress 0.1721685
## ... New best solution
## ... Procrustes: rmse 0.0009834757 max resid 0.01024635
## Run 321 stress 0.1751463
## Run 322 stress 0.1741857
## Run 323 stress 0.1747149
## Run 324 stress 0.1739038
## Run 325 stress 0.1733386
## Run 326 stress 0.1722203
## ... Procrustes: rmse 0.004457349 max resid 0.03053609
## Run 327 stress 0.175816
## Run 328 stress 0.175337
## Run 329 stress 0.1763445
## Run 330 stress 0.1726706
## Run 331 stress 0.1754963
## Run 332 stress 0.1747553
## Run 333 stress 0.173047
## Run 334 stress 0.1739535
## Run 335 stress 0.1747817
## Run 336 stress 0.1724755
## ... Procrustes: rmse 0.02017503 max resid 0.1470259
## Run 337 stress 0.1758235
## Run 338 stress 0.1735397
## Run 339 stress 0.1721884
## ... Procrustes: rmse 0.004215789 max resid 0.04444782
## Run 340 stress 0.1738999
## Run 341 stress 0.1769098
## Run 342 stress 0.1751327
## Run 343 stress 0.1730756
## Run 344 stress 0.1746818
## Run 345 stress 0.1725896
## ... Procrustes: rmse 0.02081706 max resid 0.1467707
## Run 346 stress 0.1739574
## Run 347 stress 0.1746009
## Run 348 stress 0.1737595
## Run 349 stress 0.1745198
## Run 350 stress 0.1730399
## Run 351 stress 0.1733815
## Run 352 stress 0.1735082
## Run 353 stress 0.1748932
## Run 354 stress 0.175114
## Run 355 stress 0.1733985
## Run 356 stress 0.1726582
## ... Procrustes: rmse 0.01868993 max resid 0.1450401
## Run 357 stress 0.1730416
## Run 358 stress 0.1755821
## Run 359 stress 0.176917
## Run 360 stress 0.1772515
## Run 361 stress 0.1732641
## Run 362 stress 0.1726746
## Run 363 stress 0.1750414
## Run 364 stress 0.1745283
## Run 365 stress 0.1731789
## Run 366 stress 0.1770884
## Run 367 stress 0.1742736
## Run 368 stress 0.1757257
## Run 369 stress 0.1749761
## Run 370 stress 0.1729735
## Run 371 stress 0.173791
## Run 372 stress 0.1723384
## ... Procrustes: rmse 0.0194334 max resid 0.1454182
## Run 373 stress 0.1727315
## Run 374 stress 0.1726655
## ... Procrustes: rmse 0.02124793 max resid 0.1450788
## Run 375 stress 0.1733572
## Run 376 stress 0.1723717
## ... Procrustes: rmse 0.01966895 max resid 0.144545
## Run 377 stress 0.172977
## Run 378 stress 0.1722662
## ... Procrustes: rmse 0.003896939 max resid 0.02814835
## Run 379 stress 0.1727404
## Run 380 stress 0.1731728
## Run 381 stress 0.1729195
## Run 382 stress 0.1761529
## Run 383 stress 0.1726321
## ... Procrustes: rmse 0.02534648 max resid 0.1557111
## Run 384 stress 0.174935
## Run 385 stress 0.1727799
## Run 386 stress 0.1742647
## Run 387 stress 0.1742623
## Run 388 stress 0.1723457
## ... Procrustes: rmse 0.01991658 max resid 0.1496722
## Run 389 stress 0.1722464
## ... Procrustes: rmse 0.01773663 max resid 0.1480927
## Run 390 stress 0.1761951
## Run 391 stress 0.1759265
## Run 392 stress 0.1750947
## Run 393 stress 0.1727965
## Run 394 stress 0.176099
## Run 395 stress 0.1733906
## Run 396 stress 0.1728845
## Run 397 stress 0.1725288
## ... Procrustes: rmse 0.0127141 max resid 0.1179979
## Run 398 stress 0.1752787
## Run 399 stress 0.1762001
## Run 400 stress 0.1750601
## *** No convergence -- monoMDS stopping criteria:
## 266: no. of iterations >= maxit
## 134: stress ratio > sratmaxplot(bird_nmds)
stressplot(bird_nmds)
bird species most associated with community differences
We can also look at which birds are most differentiated across communities for both the PCOA and NMDS ordinations. The below code essentially a) relates individual bird presence/absence to the final PCOA ordination, b) extracts p-values for associations with MDS1 and 2 (the first two ordination axes), and c) selects only those birds that are most differentiated at p less than or equal to 0.001. We will do so for the PCA, PCOA, and NMDS.
For all of these analyses, keep in mind that these are differentiated across plots, rather than reservations. We will explore a way to identify reservation-level differentiated species a litle bit below.
Here is the PCA:
fitpca <- envfit(p_unconc, birdcomm_bin, perm = 999)
# extract p-values for each species
fitpca_pvals <- fitpca$vectors$pvals %>%
as.data.frame() %>%
rownames_to_column("species") %>%
dplyr::rename("pvals" = ".")
# extract coordinates for species, only keep species with p-val = 0.001
fitpca_spp <- fitpca %>%
scores(., display = "vectors") %>%
as.data.frame() %>%
rownames_to_column("species") %>%
full_join(., fitpca_pvals, by = "species") %>%
filter(pvals == 0.001)47 species were identified as significantly contributing to variation between plots. These are my best attempts to break these down by typical migratory patterns for each species in Ohio (courtesy of www.allaboutbirds.org).
Migratory, breeding range in Ohio: Acadian Fly-Catcher, Baltimore Oriole, Barn Swallow, Blue-Gray Gnatcatcher, Brown Thrasher, Common Yellowthroat, Eastern Bluebird, Eastern Kingbird, Eastern Towhee, Eastern Wood-Pewee, Gray Catbird, Great-Crested Flycatcher, Hooded Warbler, House Wren, Indigo Bunting, Killdeer, Louisiana Water-Thrush, Orchard Oriole, Rose-breasted Grosbeak, Red-eyed Vireo, Scarlet Tanager, Tree Swallow, Warbling Vireo, Willow Flycatcher, Wood Thrush, Yellow Warbler, Yellow-Throated Vireo = 27 species
Year-round Ohio resident: American Crow, American Kestrel, Black-Capped Chickadee, Brown-Headed Cowbird, Cedar Waxwing, Common Grackle, Eastern Meadowlark, European Starling, Hairy Woodpecker, House Finch, House Sparrow, Mourning Dove, Northern Flicker, Pileated Woodpecker, Red-winged Blackbird, Song Sparrow, Tufted Titmouse, Turkey Vulture, White-breasted Nuthatch = 19 species
Migratory, non-breeding range in Ohio: Dark-eyed Junco = 1 species
We can repeat this procedure for the PCOA results too. Whereas PCA (above) was based on summarizing our bird community in the fewest number of axes that explain the most variance, PCOA is trying to maximize distance between plots. This is somewhat semantic in many cases, but will lead to differing representations of which birds are most associated with ordination axes.
Anyway! on to the PCOA:
# envfit() takes the output of metaMDS() and the species matrix you created
fitpcoa <- envfit(bird_pcoa, birdcomm_bin, perm = 999)
# extract p-values for each species
fitpcoa_pvals <- fitpcoa$vectors$pvals %>%
as.data.frame() %>%
rownames_to_column("species") %>%
dplyr::rename("pvals" = ".")
# extract coordinates for species, only keep species with p-val = 0.001
fit_pcoa_spp <- fitpcoa %>%
scores(., display = "vectors") %>%
as.data.frame() %>%
rownames_to_column("species") %>%
full_join(., fitpcoa_pvals, by = "species") %>%
filter(pvals == 0.001)So, now we have 40 species that contribute most to community distance between plots.
Migratory, breeding range in Ohio: Acadian Fly-Catcher, Baltimore Oriole, Barn Swallow, Blue-Gray Gnatcatcher, Brown Thrasher, Chimney Swift, Common Yellowthroat, Eastern Kingbird, Eastern Wood-Pewee, Gray Catbird, Hooded Warbler, House Wren, Indigo Bunting, Louisiana Water-Thrush, Orchard Oriole, Rose-breasted Grosbeak, Red-eyed Vireo, Scarlet Tanager, Tree Swallow, Veery, Warbling Vireo, Willow Flycatcher, Wood Thrush, Yellow Warbler, Yellow-Throated Vireo = 25 species
Year-round Ohio resident: American Robin, Brown-Headed Cowbird, Cedar Waxwing, Common Grackle, European Starling, Field Sparrow, House Finch, House Sparrow, Mourning Dove, Northern Flicker, Pileated Woodpecker, Red-winged Blackbird, Song Sparrow, White-breasted Nuthatch = 14 species
Migratory, non-breeding range in Ohio: Dark-eyed Junco = 1 species
and now the NMDS:
## species NMDS1 NMDS2 pvals
## 1 ACFL -0.6615159 -0.33993204 0.001
## 2 AMCR -0.1824643 -0.21687480 0.001
## 3 AMRO 0.1652093 0.29759255 0.001
## 4 BAOR 0.4715001 -0.14130150 0.001
## 5 BARS 0.3847375 -0.18007515 0.001
## 6 BCCH -0.1821493 -0.27828156 0.001
## 7 BGGN -0.1695426 -0.29758264 0.001
## 8 BHCO 0.2126973 -0.30699109 0.001
## 9 BRTH 0.2836090 -0.25445830 0.001
## 10 CEDW 0.4963877 -0.14194614 0.001
## 11 COGR 0.4168400 0.19899144 0.001
## 12 COYE 0.4305116 -0.45656653 0.001
## 13 DEJU -0.2774805 -0.27502180 0.001
## 14 EAKI 0.3000369 -0.15315332 0.001
## 15 EAWP -0.5238590 -0.16395515 0.001
## 16 EUST 0.5560226 0.08649803 0.001
## 17 FISP 0.1808023 -0.31446546 0.001
## 18 GRCA 0.5924405 0.01215275 0.001
## 19 HOFI 0.2812091 0.22451853 0.001
## 20 HOSP 0.2934945 0.26677074 0.001
## 21 HOWA -0.3074823 -0.54894517 0.001
## 22 HOWR 0.3033252 0.22647663 0.001
## 23 INBU 0.2723384 -0.47973161 0.001
## 24 LOWA -0.1413892 -0.35488881 0.001
## 25 MODO 0.4502863 -0.27715337 0.001
## 26 NOFL 0.3207516 -0.04470348 0.001
## 27 OROR 0.6432191 -0.04993506 0.001
## 28 PIWO -0.1199985 -0.42279736 0.001
## 29 RBGR 0.2600685 -0.34687602 0.001
## 30 REVI -0.4280763 -0.13293623 0.001
## 31 RWBL 0.6453788 -0.06557360 0.001
## 32 SCTA -0.4184239 -0.49848494 0.001
## 33 SOSP 0.6433413 -0.09766379 0.001
## 34 TRES 0.4269827 -0.11404675 0.001
## 35 WAVI 0.5914979 0.02827896 0.001
## 36 WBNU -0.3438214 -0.24357260 0.001
## 37 WIFL 0.3800642 -0.27789620 0.001
## 38 WOTH -0.1567631 -0.52087264 0.001
## 39 YEWA 0.6798322 -0.20922757 0.001
## 40 YTVI -0.1698554 -0.34639158 0.001So, now we have 41 species that contribute most to community distance between plots.
Migratory, breeding range in Ohio: Acadian Fly-Catcher, Baltimore Oriole, Barn Swallow, Blue-Gray Gnatcatcher, Brown Thrasher, Common Yellowthroat, Eastern Kingbird, Eastern Wood-Pewee, Gray Catbird, Hooded Warbler, House Wren, Indigo Bunting, Louisiana Water-Thrush, Orchard Oriole, Rose-breasted Grosbeak, Red-eyed Vireo, Scarlet Tanager, Tree Swallow, Veery, Warbling Vireo, Willow Flycatcher, Wood Thrush, Yellow Warbler, Yellow-Throated Vireo
Year-round Ohio resident: American Robin, Black-Capped Chickadee, Brown-Headed Cowbird, Canada Goose, Cedar Waxwing, Common Grackle, European Starling, Field Sparrow, House Finch, House Sparrow, Mourning Dove, Northern Flicker, Pileated Woodpecker, Red-winged Blackbird, Song Sparrow, White-breasted Nuthatch
Migratory, non-breeding range in Ohio: Dark-eyed Junco
While these three sets of species are pretty similar, we can also take a look at which species are different across approaches.
outersect(fit_spp$species,fit_pcoa_spp$species,fitpca_spp$species)## [1] "CHSW" "EATO" "VEER" "AMGO" "AMKE" "CANG" "EABL" "EAME" "GCFL" "KILL" "TUVU"So American Crow, American Kestrel, Canada Goose, Chimney Swift, Eastern Bluebird, Eastern Meadowlark, Easter Towhee, Great-Crested Flycatcher, Hairy Woodpecker, Killdeer, Tufted Titmouse, and Turkey Vulture. At least to me, nothing particularly sticks out as a commonality across these species, but it may make just as much sense to simply pool these species across approaches:
sort(unique(c(fit_spp$species,fit_pcoa_spp$species,fitpca_spp$species)))## [1] "ACFL" "AMCR" "AMGO" "AMKE" "AMRO" "BAOR" "BARS" "BCCH" "BGGN" "BHCO" "BRTH" "CANG" "CEDW" "CHSW" "COGR" "COYE" "DEJU" "EABL" "EAKI" "EAME" "EATO" "EAWP" "EUST" "FISP" "GCFL" "GRCA" "HAWO" "HOFI" "HOSP" "HOWA" "HOWR" "INBU" "KILL" "LOWA" "MODO" "NOFL" "NRWS" "OROR" "PIWO" "RBGR" "REVI" "RWBL"
## [43] "SCTA" "SOSP" "TRES" "TUVU" "VEER" "WAVI" "WBNU" "WIFL" "WOTH" "YEWA" "YTVI"This gives us a list of 52 bird species that seem to be particularly differentiated/variable across plots. That is precisely half of the total number of species in the overall dataset.
As we proceed, we will begin layering on our environmental datasets using three main approaches: constrained ordination, linear models of species richness, and multi-species occupancy models.
constrained
Whereas the ordination approaches above were focused on unconstrained ordination - in this case, exploring how species are arranged between sample sites - there is an entire suite of approaches that address how environmental gradients shape community assemblages among sample sites. These are called constrained ordination approaches. I am personally a huge proponent of redundancy analysis (RDA) as a flexible constrained ordination approach. Taken from the GUSTA ME tutorial on RDA:
“RDA can also be considered a constrained version of principal components analysis (PCA), wherein canonical axes - built from linear combinations of response variables - must also be linear combinations of the explanatory variables (i.e. fitted by [multiple linear regression]). The RDA approach generates one ordination in the space defined by the matrix of response variables and another in the space defined by the matrix of explanatory variables.”
So, in essence, the ordination plot produced will depict an optimal arrangement of environmental axes (RDA1, 2, etc.) that maximally explains variation in the species dataset.
We are going to use stepwise selection to explore which PCAP covariates most explain plot-level differences.
p<-rda(birdcomm_bin~.,vegecomm_drop,scale=T)
fit<-envfit(p,vegecomm_drop,scale=T)
p_0<-rda(birdcomm_bin~1,vegecomm_drop,scale=T)
p_1<-rda(birdcomm_bin~.,vegecomm_drop,scale=T)
set.seed(12345)
p_step<- ordistep(p_0, scope = formula(p_1))##
## Start: birdcomm_bin ~ 1
##
## Df AIC F Pr(>F)
## + tol_rel 1 756.52 6.6214 0.005 **
## + svp 1 759.73 3.3565 0.005 **
## + sens_rel 1 760.09 2.9905 0.005 **
## + bryo 1 760.14 2.9425 0.005 **
## + shrub 1 760.70 2.3844 0.005 **
## + dicot 1 760.72 2.3572 0.005 **
## + sens_rel_notree 1 760.89 2.1925 0.005 **
## + cyper 1 761.26 1.8221 0.005 **
## + canopy 1 761.41 1.6719 0.005 **
## + subcanopy 1 761.12 1.9574 0.010 **
## + inv 1 760.60 2.4826 0.060 .
## + hydro 1 761.74 1.3399 0.100 .
## + small_tree 1 761.87 1.2124 0.120
## + hydro_rel 1 762.00 1.0866 0.290
## + ap 1 762.19 0.8948 0.500
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Step: birdcomm_bin ~ tol_rel
##
## Df AIC F Pr(>F)
## - tol_rel 1 761.09 6.6214 0.005 **
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Df AIC F Pr(>F)
## + dicot 1 756.28 2.2181 0.005 **
## + svp 1 756.37 2.1250 0.005 **
## + cyper 1 756.72 1.7772 0.005 **
## + shrub 1 757.08 1.4239 0.040 *
## + hydro 1 757.11 1.3912 0.050 *
## + canopy 1 757.22 1.2845 0.055 .
## + bryo 1 757.02 1.4847 0.070 .
## + inv 1 756.58 1.9146 0.080 .
## + small_tree 1 757.29 1.2147 0.130
## + subcanopy 1 757.41 1.0943 0.310
## + hydro_rel 1 757.47 1.0400 0.345
## + sens_rel_notree 1 757.46 1.0427 0.395
## + sens_rel 1 757.50 1.0083 0.435
## + ap 1 757.66 0.8463 0.600
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Step: birdcomm_bin ~ tol_rel + dicot
##
## Df AIC F Pr(>F)
## - dicot 1 756.52 2.2181 0.005 **
## - tol_rel 1 760.72 6.4527 0.005 **
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Df AIC F Pr(>F)
## + svp 1 756.32 1.9190 0.005 **
## + cyper 1 756.36 1.8840 0.005 **
## + bryo 1 756.74 1.5059 0.060 .
## + inv 1 756.33 1.9124 0.075 .
## + canopy 1 756.95 1.3002 0.080 .
## + small_tree 1 757.04 1.2114 0.140
## + shrub 1 757.17 1.0897 0.285
## + sens_rel_notree 1 757.22 1.0397 0.340
## + subcanopy 1 757.23 1.0279 0.390
## + hydro_rel 1 757.25 1.0092 0.395
## + sens_rel 1 757.23 1.0254 0.410
## + ap 1 757.40 0.8590 0.550
## + hydro 1 757.35 0.9127 0.695
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Step: birdcomm_bin ~ tol_rel + dicot + svp
##
## Df AIC F Pr(>F)
## - dicot 1 756.37 2.0114 0.010 **
## - svp 1 756.28 1.9190 0.005 **
## - tol_rel 1 759.74 5.3757 0.005 **
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Df AIC F Pr(>F)
## + cyper 1 756.50 1.7813 0.005 **
## + canopy 1 756.97 1.3150 0.070 .
## + inv 1 756.35 1.9287 0.075 .
## + bryo 1 756.98 1.3094 0.165
## + small_tree 1 757.12 1.1738 0.200
## + hydro_rel 1 757.27 1.0211 0.340
## + shrub 1 757.29 1.0049 0.405
## + subcanopy 1 757.29 1.0052 0.425
## + sens_rel_notree 1 757.30 0.9966 0.440
## + ap 1 757.43 0.8715 0.535
## + sens_rel 1 757.34 0.9545 0.605
## + hydro 1 757.38 0.9145 0.610
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Step: birdcomm_bin ~ tol_rel + dicot + svp + cyper
##
## Df AIC F Pr(>F)
## - cyper 1 756.32 1.7813 0.010 **
## - svp 1 756.36 1.8160 0.005 **
## - dicot 1 756.83 2.2826 0.005 **
## - tol_rel 1 759.99 5.4185 0.005 **
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Df AIC F Pr(>F)
## + inv 1 756.50 1.9396 0.080 .
## + bryo 1 757.12 1.3309 0.085 .
## + canopy 1 757.18 1.2729 0.100 .
## + small_tree 1 757.31 1.1473 0.260
## + shrub 1 757.41 1.0509 0.320
## + hydro_rel 1 757.41 1.0481 0.385
## + sens_rel_notree 1 757.45 1.0108 0.400
## + subcanopy 1 757.47 0.9963 0.525
## + sens_rel 1 757.51 0.9578 0.560
## + hydro 1 757.54 0.9208 0.600
## + ap 1 757.70 0.7733 0.690
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1p_step## Call: rda(formula = birdcomm_bin ~ tol_rel + dicot + svp + cyper, data = vegecomm_drop, scale = T)
##
## Inertia Proportion Rank
## Total 103.00000 1.00000
## Constrained 7.61479 0.07393 4
## Unconstrained 95.38521 0.92607 102
## Inertia is correlations
##
## Eigenvalues for constrained axes:
## RDA1 RDA2 RDA3 RDA4
## 4.483 1.406 0.991 0.735
##
## Eigenvalues for unconstrained axes:
## PC1 PC2 PC3 PC4 PC5 PC6 PC7 PC8
## 6.458 3.725 3.535 3.198 3.030 2.685 2.578 2.387
## (Showing 8 of 102 unconstrained eigenvalues)anova(p_step)## Permutation test for rda under reduced model
## Permutation: free
## Number of permutations: 999
##
## Model: rda(formula = birdcomm_bin ~ tol_rel + dicot + svp + cyper, data = vegecomm_drop, scale = T)
## Df Variance F Pr(>F)
## Model 4 7.615 3.1733 0.001 ***
## Residual 159 95.385
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1anova(p_step,by="term")## Permutation test for rda under reduced model
## Terms added sequentially (first to last)
## Permutation: free
## Number of permutations: 999
##
## Model: rda(formula = birdcomm_bin ~ tol_rel + dicot + svp + cyper, data = vegecomm_drop, scale = T)
## Df Variance F Pr(>F)
## tol_rel 1 4.045 6.7420 0.001 ***
## dicot 1 1.345 2.2417 0.001 ***
## svp 1 1.157 1.9283 0.003 **
## cyper 1 1.069 1.7813 0.006 **
## Residual 159 95.385
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1anova(p_step,by="axis")## Permutation test for rda under reduced model
## Forward tests for axes
## Permutation: free
## Number of permutations: 999
##
## Model: rda(formula = birdcomm_bin ~ tol_rel + dicot + svp + cyper, data = vegecomm_drop, scale = T)
## Df Variance F Pr(>F)
## RDA1 1 4.483 7.4724 0.001 ***
## RDA2 1 1.406 2.3436 0.001 ***
## RDA3 1 0.991 1.6527 0.007 **
## RDA4 1 0.735 1.2246 0.095 .
## Residual 159 95.385
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1ordiplot(p_step,scaling=3)
So this tells us that the most explanatory group of PCAP covariates is relative cover of tolerant plants (tol_rel), cover of dicots, cover of seedless vascular plants (svp), and cover of sedges (cyper). Rather conveniently, the relative cover of dicots, seedless vascular plants, and sedges is negatively associated with the relative cover of tolerant plant species. As a quick caveat, these four variables only explain about 7.3% of the variation, as seen here (red is our environmental axes, black is our unconstrained axes):
constrained_eig <- p_step$CCA$eig/p_step$tot.chi*100
unconstrained_eig <- p_step$CA$eig/p_step$tot.chi*100
expl_var <- c(constrained_eig, unconstrained_eig)
barplot (expl_var[1:20], col = c(rep ('red', length (constrained_eig)), rep ('black', length (unconstrained_eig))),
las = 2, ylab = '% variation')
This is not that terribly uncommon for redundancy analyses - It is often difficult to find environmental gradients that explain all of the variance in large datasets like ours. Biologically, we likely do not necessarily expect that a single environmental gradient alone drives all of the variation in bird presence/absence, and intuitively know that we also have to worry about observer effects, weather, and other unmeasured characteristics. Regardless of any of the above, the ordination is significant - as are the relationships with these covariates, so we can safely say that nearly 10% of the variation in our bird species dataset is explainable by these 4 PCAP covariates!
We can also explore bird species most significantly associated with our constrained ordination axes, too.
# envfit() takes the output of metaMDS() and the species matrix you created
fitconcrda <- envfit(p_step, birdcomm_bin, perm = 999)
# extract p-values for each species
fitconcrda_pvals <- fitconcrda$vectors$pvals %>%
as.data.frame() %>%
rownames_to_column("species") %>%
dplyr::rename("pvals" = ".")
# extract coordinates for species, only keep species with p-val = 0.001
fitconcrda_spp <- fitconcrda %>%
scores(., display = "vectors") %>%
as.data.frame() %>%
rownames_to_column("species") %>%
full_join(., fitconcrda_pvals, by = "species") %>%
filter(pvals == 0.001)OK - this is getting towards our main analytical goals! Let’s sort by scores for that first environmental axis and plot these in order:
fitconcrda_spp_order<-fitconcrda_spp[order(fitconcrda_spp$RDA1),]
fitconcrda_spp_order$species<-factor(fitconcrda_spp_order$species,levels=fitconcrda_spp_order$species)
ggplot(aes(x=RDA1,y=species),data=fitconcrda_spp_order)+geom_point()+theme_classic(base_size=16)
This is a pretty neat plot - Acadian Flycatcher, Scarlet Tanager, Eastern Wood-Pewee, White-breasted Nuthatch, Hooded Warbler, Red-Eyed Vireo, Dark-Eyed Junco, Wood Thrush, Hairy Woodpeckers, and Lousiana Water-Thrush all appear to be associated with higher sedge, seedless vascular plant, and dicot cover and lower cover of tolerant plants. In contrast, Yellow Warblers, Orchard Orioles, Red-Winged Blackbirds, European Starlings, Warbling Vireos, Song Sparrows, Grey Catbirds, Tree Swallows, Mourning Doves, and Cedar Waxwings are all associated with more tolerant plants.
constrained - sites
While we have already explored environmental associations, we might also want to know what species are most differentiated between reservations. We will then use similar code to that used above to explore species that are most differentiated between reservations.
p_reservations<-rda(birdcomm_bin~reservatio,CMP_bbs_pcap_join_pcaponly_pcapsum,scale=T)
ordiplot(p_reservations,scaling=3)
reservations_fitbirds<-envfit(p_reservations,birdcomm_bin)
reservations_fit<-envfit(p_reservations~reservatio,data=CMP_bbs_pcap_join_pcaponly_pcapsum,scale=T)
CMP_bbs_pcap_join_pcaponly_pcapsum$site<-rownames(CMP_bbs_pcap_join_pcaponly_pcapsum)
plot_df_reservations <- scores(p_reservations, display = "sites") %>%
as.data.frame() %>%
rownames_to_column("site") %>%
full_join(CMP_bbs_pcap_join_pcaponly_pcapsum, by = "site")
plot_rda_reservations <- ggplot(plot_df_reservations, aes(x = RDA1, y = RDA2, color = reservatio)) +
geom_point(size = 3, alpha = 0.8) +
stat_ellipse(linetype = 2, size = 1) +
labs(title = "RDA")
plot_rda_reservations## Too few points to calculate an ellipse
## Too few points to calculate an ellipse
## Too few points to calculate an ellipse
## Too few points to calculate an ellipse
## Too few points to calculate an ellipse
## Too few points to calculate an ellipse
## Too few points to calculate an ellipse## Warning: Removed 7 row(s) containing missing values (geom_path).
fitreservations_pvals <- reservations_fitbirds$vectors$pvals %>%
as.data.frame() %>%
rownames_to_column("species") %>%
dplyr::rename("pvals" = ".")
# extract coordinates for species, only keep species with p-val = 0.001
fit_reservations_spp <- reservations_fitbirds %>%
scores(., display = "vectors") %>%
as.data.frame() %>%
rownames_to_column("species") %>%
full_join(., fitreservations_pvals, by = "species") %>%
filter(pvals == 0.001)
fit_reservations_spp## species RDA1 RDA2 pvals
## 1 ACFL -0.67536570 0.09819212 0.001
## 2 AMCR -0.37692846 -0.20456090 0.001
## 3 AMRE -0.07400135 -0.35850489 0.001
## 4 AMRO 0.26495737 0.23369097 0.001
## 5 BAOR 0.41712060 0.17739155 0.001
## 6 BARS 0.36240761 -0.30012415 0.001
## 7 BCCH -0.24236822 0.20415669 0.001
## 8 BRCR -0.28606094 -0.09945096 0.001
## 9 BRTH 0.25037094 -0.24226047 0.001
## 10 CEDW 0.41740795 -0.28403651 0.001
## 11 COGR 0.51954277 0.18156226 0.001
## 12 COYE 0.09098344 -0.54427390 0.001
## 13 DEJU -0.37988349 -0.03903739 0.001
## 14 EATO -0.12707725 -0.27487048 0.001
## 15 EAWP -0.44389257 0.47561627 0.001
## 16 EUST 0.55101057 -0.40645368 0.001
## 17 FISP -0.02103263 -0.34546312 0.001
## 18 GRCA 0.46149939 -0.09909873 0.001
## 19 HAWO -0.14783078 0.26640131 0.001
## 20 HOFI 0.39273807 -0.02749710 0.001
## 21 HOSP 0.47906450 0.11992733 0.001
## 22 HOWA -0.56788748 -0.25735546 0.001
## 23 HOWR 0.32157965 0.13443633 0.001
## 24 INBU 0.14260564 -0.26880856 0.001
## 25 LOWA -0.29872249 -0.25058885 0.001
## 26 MODO 0.36998816 -0.15778881 0.001
## 27 NOFL 0.30988119 0.16027358 0.001
## 28 NRWS 0.44336152 0.10254562 0.001
## 29 OROR 0.54726150 -0.40264845 0.001
## 30 OVEN -0.27119395 -0.22211081 0.001
## 31 RBGU 0.36218342 -0.09026430 0.001
## 32 REVI -0.29178905 0.32096716 0.001
## 33 RWBL 0.50726734 -0.28022387 0.001
## 34 SCTA -0.54059278 0.06735648 0.001
## 35 SOSP 0.46731737 -0.16498066 0.001
## 36 TRES 0.44251225 -0.11105510 0.001
## 37 VEER -0.25041023 -0.17876454 0.001
## 38 WAVI 0.54355113 -0.29016743 0.001
## 39 WBNU -0.45277798 0.32182366 0.001
## 40 WIFL 0.28267826 -0.33559442 0.001
## 41 WOTH -0.35735577 -0.01891986 0.001
## 42 YEWA 0.50674786 -0.44480078 0.001
## 43 YTVI -0.35813458 -0.21812489 0.001fit_reservations <- reservations_fit %>%
scores(., display = "factors") %>%
as.data.frame() %>%
rownames_to_column("reservations")
rda_reservation_plot_new <- ggplot(plot_df_reservations, aes(x = RDA1, y = RDA2)) +
coord_fixed() +
geom_point(aes(color = reservatio), size = 3, alpha = 0.8) +
geom_segment(data = fit_reservations_spp, aes(x = 0, xend = RDA1, y = 0, yend = RDA2),
col = "black") +
geom_text(data = fit_reservations_spp, aes(label = species))+
geom_text(data=as.data.frame(fit_reservations),aes(label=reservations))
rda_reservation_plot_new
Interesting! So when we consider reservation-of-origin, Hinckley seems to have a pretty unique bird community overall, and groups negatively along RDA1. In contrast, Brookside, Ohio and Erie Canal, and Washington all group positively along RDA1 - and are all quite heavily urbanized.
While I personally think the environmental gradient analyses are a bit more interesting, this does highlight reservations that are more or less distinct in bird community space. We can also plot just the reservations alone too:
ggplot(plot_df_reservations, aes(x = RDA1, y = RDA2)) +
coord_fixed() +
geom_text(data=as.data.frame(fit_reservations),aes(label=reservations))+xlim(-2,4)+ylim(-4,2)
Now pRDA - with variance partitioning.
As a final sanity check, it can sometimes be difficult to disentangle geographic effects from environmental effects in these sorts of analyses. For example, what if (hypothetically) cover of seedless vascular plants is negatively associated with geographic coordinates? Thus, we might want to briefly check if the variance explained by our environmental axes is (or is not) confounded with geography. To do so, we will use the varpart() function, which partitions variance between model components:
vegecomm_latlong<-CMP_bbs_pcap_join_pcaponly_pcapsum[,2:26]
pcap_geog<-varpart(birdcomm_bin,~ tol_rel + svp + dicot + cyper,~longitude+latitude,data=vegecomm_latlong)
plot(pcap_geog)
So - most of this variation is associated with these top covariates, and NOT with latitude and longitude alone. Good to see - this means we can be confident that these are actual environmental effects!
community richness metrics
Next, we can explore the relationships between species richness and plot characteristics. We will start with plot-level explorations.
plot-level analyses
plot-level richness and various indicators
Now, we will produce some exploratory models and plots of richness against basic indicators of vegetative quality, including the Floristic Quality Indicator (FQI) and the Vegetative Index of Biotic Integrity (VIBI).
S <- specnumber(birdcomm_bin)
vegecomm_richness<-cbind(vegecomm,S)
cor(vegecomm_richness)## vibifq vibifqnotrees carex cyper dicot shade shrub hydro svp ap fqi bryo hydro_rel hydro_rel_notrees sens_rel sens_rel_notree tol_rel tol_rel_notree inv small_tree subcanopy
## vibifq 1.000000000 0.903466506 0.302129838 0.268662618 0.44969374 0.643467211 0.42091385 0.094920118 0.52435545 -0.201584848 9.191368e-01 0.227249850 0.14141059 0.232333080 0.615943742 0.517291454 -0.82095096 -0.815660534 -0.204274958 -8.896552e-03 0.352987107
## vibifqnotrees 0.903466506 1.000000000 0.315594820 0.286141876 0.57241339 0.762505137 0.44919931 0.183868172 0.57125442 -0.163682237 9.157120e-01 0.179820785 0.23974361 0.250108233 0.401409023 0.573829627 -0.66469643 -0.758049585 -0.150119819 4.913246e-03 0.428733969
## carex 0.302129838 0.315594820 1.000000000 0.990155702 0.46569343 0.543313784 0.14791367 0.173971056 0.32038084 -0.207379567 4.247642e-01 -0.088332875 -0.01680755 -0.001503659 -0.010702821 0.022201528 -0.08217802 -0.040658648 -0.068732762 4.256074e-02 0.098208268
## cyper 0.268662618 0.286141876 0.990155702 1.000000000 0.45810906 0.519855636 0.12970354 0.173572197 0.30826804 -0.214455749 3.923608e-01 -0.092992405 -0.01088781 0.001649045 -0.026296645 0.010052307 -0.04666473 -0.008993867 -0.072806020 4.618045e-02 0.105018429
## dicot 0.449693738 0.572413388 0.465693427 0.458109060 1.00000000 0.900954560 0.38443818 0.501981435 0.48161327 -0.120998160 6.868442e-01 -0.078551692 0.22182199 0.148547049 -0.038371772 0.153124228 -0.09258616 -0.111832162 -0.074575448 -1.060284e-01 0.325216371
## shade 0.643467211 0.762505137 0.543313784 0.519855636 0.90095456 1.000000000 0.40882995 0.378278862 0.58861737 -0.156097164 8.368920e-01 0.020607143 0.20374398 0.157734262 0.086603981 0.266173043 -0.31108544 -0.344672034 -0.059325341 8.208519e-03 0.371162537
## shrub 0.420913850 0.449199309 0.147913674 0.129703536 0.38443818 0.408829955 1.00000000 0.588835377 0.29470052 -0.101698461 4.747080e-01 0.085207890 0.49360022 0.536538046 0.162917290 0.160089484 -0.22690227 -0.273227137 -0.107336552 -8.936686e-02 0.255636573
## hydro 0.094920118 0.183868172 0.173971056 0.173572197 0.50198144 0.378278862 0.58883538 1.000000000 0.17700645 -0.039802022 2.285981e-01 -0.081936068 0.54679684 0.488784912 -0.086929406 -0.061281690 0.10928036 0.079007210 0.003619817 -3.202940e-02 0.204315768
## svp 0.524355447 0.571254416 0.320380841 0.308268040 0.48161327 0.588617368 0.29470052 0.177006451 1.00000000 -0.147271791 5.999081e-01 0.244051454 0.15150342 0.153627835 0.190425176 0.258147075 -0.27430569 -0.331328168 -0.113257408 -9.341899e-02 0.268499329
## ap -0.201584848 -0.163682237 -0.207379567 -0.214455749 -0.12099816 -0.156097164 -0.10169846 -0.039802022 -0.14727179 1.000000000 -2.573751e-01 0.006014263 -0.03413775 -0.021095091 -0.055016524 -0.012339374 0.12954362 0.125069560 0.033126004 -7.653394e-02 -0.095063187
## fqi 0.919136798 0.915712021 0.424764247 0.392360839 0.68684421 0.836892017 0.47470799 0.228598128 0.59990814 -0.257375101 1.000000e+00 0.164507181 0.15874283 0.185026315 0.395205020 0.452238507 -0.64079294 -0.658950801 -0.158986465 -1.659484e-05 0.356068079
## bryo 0.227249850 0.179820785 -0.088332875 -0.092992405 -0.07855169 0.020607143 0.08520789 -0.081936068 0.24405145 0.006014263 1.645072e-01 1.000000000 0.02837391 0.046954292 0.289335075 0.286497590 -0.32259935 -0.349022222 -0.101490062 3.167128e-01 -0.004418917
## hydro_rel 0.141410593 0.239743608 -0.016807547 -0.010887815 0.22182199 0.203743981 0.49360022 0.546796836 0.15150342 -0.034137749 1.587428e-01 0.028373905 1.00000000 0.899920418 -0.091789549 -0.064596088 -0.12854831 -0.192750547 -0.048956575 6.854047e-02 0.412103596
## hydro_rel_notrees 0.232333080 0.250108233 -0.001503659 0.001649045 0.14854705 0.157734262 0.53653805 0.488784912 0.15362783 -0.021095091 1.850263e-01 0.046954292 0.89992042 1.000000000 0.061436676 -0.073288844 -0.24903476 -0.272482935 -0.076359781 2.939615e-02 0.387893166
## sens_rel 0.615943742 0.401409023 -0.010702821 -0.026296645 -0.03837177 0.086603981 0.16291729 -0.086929406 0.19042518 -0.055016524 3.952050e-01 0.289335075 -0.09178955 0.061436676 1.000000000 0.529677884 -0.55933338 -0.537226545 -0.075091400 -7.745522e-03 0.003765294
## sens_rel_notree 0.517291454 0.573829627 0.022201528 0.010052307 0.15312423 0.266173043 0.16008948 -0.061281690 0.25814708 -0.012339374 4.522385e-01 0.286497590 -0.06459609 -0.073288844 0.529677884 1.000000000 -0.38059052 -0.445572344 0.002392817 2.451389e-02 0.060370376
## tol_rel -0.820950963 -0.664696430 -0.082178015 -0.046664726 -0.09258616 -0.311085443 -0.22690227 0.109280360 -0.27430569 0.129543623 -6.407929e-01 -0.322599345 -0.12854831 -0.249034760 -0.559333375 -0.380590517 1.00000000 0.923531143 0.266063420 -1.012597e-01 -0.288508481
## tol_rel_notree -0.815660534 -0.758049585 -0.040658648 -0.008993867 -0.11183216 -0.344672034 -0.27322714 0.079007210 -0.33132817 0.125069560 -6.589508e-01 -0.349022222 -0.19275055 -0.272482935 -0.537226545 -0.445572344 0.92353114 1.000000000 0.214269092 -9.522446e-02 -0.323032669
## inv -0.204274958 -0.150119819 -0.068732762 -0.072806020 -0.07457545 -0.059325341 -0.10733655 0.003619817 -0.11325741 0.033126004 -1.589865e-01 -0.101490062 -0.04895657 -0.076359781 -0.075091400 0.002392817 0.26606342 0.214269092 1.000000000 6.041087e-02 -0.073180439
## small_tree -0.008896552 0.004913246 0.042560741 0.046180446 -0.10602839 0.008208519 -0.08936686 -0.032029403 -0.09341899 -0.076533935 -1.659484e-05 0.316712780 0.06854047 0.029396146 -0.007745522 0.024513894 -0.10125967 -0.095224459 0.060410871 1.000000e+00 -0.150833647
## subcanopy 0.352987107 0.428733969 0.098208268 0.105018429 0.32521637 0.371162537 0.25563657 0.204315768 0.26849933 -0.095063187 3.560681e-01 -0.004418917 0.41210360 0.387893166 0.003765294 0.060370376 -0.28850848 -0.323032669 -0.073180439 -1.508336e-01 1.000000000
## canopy 0.203002590 0.102153956 -0.102731833 -0.111819839 -0.02406017 0.020867114 0.10411367 0.040288929 0.07700871 -0.055629052 1.434842e-01 0.104404368 0.04753125 0.107445181 0.264437686 0.096605726 -0.30995440 -0.253543359 -0.102465127 5.548526e-02 0.015482845
## S -0.180202756 -0.094422012 -0.123628799 -0.116601889 0.14178910 0.010138092 -0.09358624 0.123776993 -0.11846054 0.087142662 -1.110069e-01 -0.190998403 -0.03126316 -0.103537313 -0.096502774 -0.027700516 0.23555710 0.198107496 0.047266522 2.431789e-02 -0.094543020
## canopy S
## vibifq 0.20300259 -0.18020276
## vibifqnotrees 0.10215396 -0.09442201
## carex -0.10273183 -0.12362880
## cyper -0.11181984 -0.11660189
## dicot -0.02406017 0.14178910
## shade 0.02086711 0.01013809
## shrub 0.10411367 -0.09358624
## hydro 0.04028893 0.12377699
## svp 0.07700871 -0.11846054
## ap -0.05562905 0.08714266
## fqi 0.14348423 -0.11100691
## bryo 0.10440437 -0.19099840
## hydro_rel 0.04753125 -0.03126316
## hydro_rel_notrees 0.10744518 -0.10353731
## sens_rel 0.26443769 -0.09650277
## sens_rel_notree 0.09660573 -0.02770052
## tol_rel -0.30995440 0.23555710
## tol_rel_notree -0.25354336 0.19810750
## inv -0.10246513 0.04726652
## small_tree 0.05548526 0.02431789
## subcanopy 0.01548284 -0.09454302
## canopy 1.00000000 -0.08082640
## S -0.08082640 1.00000000S_vibi<-lm(S~vibifq,vegecomm_richness)
summary(S_vibi)##
## Call:
## lm(formula = S ~ vibifq, data = vegecomm_richness)
##
## Residuals:
## Min 1Q Median 3Q Max
## -17.0444 -4.2186 -0.2331 3.8289 19.0802
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 28.41399 1.07526 26.425 <2e-16 ***
## vibifq -0.04908 0.02105 -2.332 0.0209 *
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 6.217 on 162 degrees of freedom
## Multiple R-squared: 0.03247, Adjusted R-squared: 0.0265
## F-statistic: 5.437 on 1 and 162 DF, p-value: 0.02094S_fqi<-lm(S~fqi,vegecomm_richness)
summary(S_fqi)##
## Call:
## lm(formula = S ~ fqi, data = vegecomm_richness)
##
## Residuals:
## Min 1Q Median 3Q Max
## -17.1875 -4.4209 -0.1052 3.9411 19.4469
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 28.1082 1.4444 19.461 <2e-16 ***
## fqi -0.1079 0.0759 -1.422 0.157
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 6.282 on 162 degrees of freedom
## Multiple R-squared: 0.01232, Adjusted R-squared: 0.006226
## F-statistic: 2.021 on 1 and 162 DF, p-value: 0.157ggplot(aes(x=vibifq,y=S),data=vegecomm_richness)+geom_point()+geom_smooth()
ggplot(aes(x=fqi,y=S),data=vegecomm_richness)+geom_point()+geom_smooth()
Somewhat unexpectedly, this seems to suggest that higher VIBI scores actually produce slightly less species-rich communities. What could be causing this? Well, to start - perhaps the number of surveys for each plot may be influencing this.
We can generate a summary of the number of surveys for each plot.
pcap_surveynumber<-CMP_bbs_pcap_join_pcaponly %>%
group_by(simp_plot) %>%
summarise(n=n_distinct(date)
)
vegecomm_richness$surveycount<-pcap_surveynumber$n
specaccum(birdcomm_bin,w=pcap_surveynumber$n,method="random")## Species Accumulation Curve
## Accumulation method: random, with 100 permutations, weighted
## Call: specaccum(comm = birdcomm_bin, method = "random", w = pcap_surveynumber$n)
##
##
## Sites 1.0000 2.0000 3.0000 4.0000 5.0000 6.0000 7.0000 8.0000 9.0000 10.0000 11.000 12.0000 13.0000 14.0000 15.0000 16.0000 17.0000 18.0000 19.0000 20.0000 21.0000 22.0000 23.0000 24.0000 25.0000 26.0000 27.0000 28.0000 29.0000 30.0000 31.0000 32.0000 33.0000 34.0000
## Effort 4.6646 9.3293 13.9939 18.6585 23.3232 27.9878 32.6524 37.3171 41.9817 46.6463 51.311 55.9756 60.6402 65.3049 69.9695 74.6341 79.2988 83.9634 88.6280 93.2927 97.9573 102.6220 107.2866 111.9512 116.6159 121.2805 125.9451 130.6098 135.2744 139.9390 144.6037 149.2683 153.9329 158.5976
## Richness NA 37.0601 43.9160 48.4527 52.0466 55.2858 58.1210 60.0900 61.7656 63.1583 64.629 66.1469 67.4490 68.7554 69.9239 70.9281 71.7678 72.6455 73.4860 74.3025 74.8494 75.6023 76.3352 76.8898 77.5041 77.9124 78.4497 79.0296 79.4300 79.9113 80.3283 80.6705 81.0429 81.4608
## sd NA 5.4585 5.3625 5.5323 5.7458 5.6377 5.4931 5.4603 5.1609 5.0987 4.744 4.5960 4.5160 4.2286 4.2348 4.1589 4.0009 4.1083 4.1018 3.9875 4.0215 3.9500 3.7498 3.6637 3.6386 3.5648 3.5230 3.3841 3.3204 3.2723 3.2496 3.2391 3.1904 3.1615
##
## Sites 35.0000 36.0000 37.0000 38.0000 39.0000 40.0000 41.0000 42.0000 43.0000 44.0000 45.0000 46.0000 47.0000 48.0000 49.0000 50.0000 51.0000 52.0000 53.0000 54.0000 55.0000 56.0000 57.0000 58.0000 59.0000 60.0000 61.0000 62.0000 63.0000 64.0000 65.0000 66.0000
## Effort 163.2622 167.9268 172.5915 177.2561 181.9207 186.5854 191.2500 195.9146 200.5793 205.2439 209.9085 214.5732 219.2378 223.9024 228.5671 233.2317 237.8963 242.5610 247.2256 251.8902 256.5549 261.2195 265.8841 270.5488 275.2134 279.8780 284.5427 289.2073 293.8720 298.5366 303.2012 307.8659
## Richness 81.9161 82.3425 82.7161 83.0964 83.4936 83.8933 84.2293 84.5283 84.8701 85.2724 85.5387 85.8263 86.1494 86.4410 86.7592 87.0225 87.3165 87.6131 87.8547 88.0760 88.3363 88.5027 88.6978 88.9875 89.2456 89.4841 89.7261 89.8954 90.0932 90.3161 90.4846 90.6542
## sd 3.1399 3.0622 3.0889 2.9418 2.8768 2.8643 2.8771 2.9532 2.9370 2.9479 2.9663 2.9171 2.8620 2.8931 2.8170 2.7562 2.7547 2.8049 2.7477 2.7540 2.7294 2.7151 2.7060 2.7082 2.6848 2.6323 2.6701 2.6262 2.6682 2.6914 2.6771 2.6663
##
## Sites 67.0000 68.0000 69.0000 70.0000 71.0000 72.0000 73.0000 74.0000 75.0000 76.0000 77.0000 78.0000 79.0000 80.0000 81.0000 82.0000 83.0000 84.0000 85.0000 86.0000 87.0000 88.0000 89.0000 90.0000 91.0000 92.0000 93.0000 94.0000 95.0000 96.0000 97.0000 98.0000
## Effort 312.5305 317.1951 321.8598 326.5244 331.1890 335.8537 340.5183 345.1829 349.8476 354.5122 359.1768 363.8415 368.5061 373.1707 377.8354 382.5000 387.1646 391.8293 396.4939 401.1585 405.8232 410.4878 415.1524 419.8171 424.4817 429.1463 433.8110 438.4756 443.1402 447.8049 452.4695 457.1341
## Richness 90.8779 91.0779 91.2097 91.3421 91.5238 91.6833 91.8987 92.1100 92.2622 92.3996 92.4958 92.6370 92.8291 92.9901 93.1603 93.3885 93.5503 93.6618 93.8185 93.9687 94.1282 94.3118 94.4901 94.6310 94.8387 94.9822 95.1005 95.2411 95.3443 95.4566 95.5830 95.7483
## sd 2.6270 2.6534 2.6440 2.5938 2.6117 2.6135 2.6498 2.6679 2.6908 2.6421 2.6006 2.5911 2.6056 2.6737 2.6933 2.5987 2.5340 2.5391 2.5614 2.5676 2.5888 2.5099 2.5361 2.4771 2.4544 2.4693 2.4538 2.4307 2.4065 2.3788 2.3770 2.4033
##
## Sites 99.0000 100.0000 101.0000 102.0000 103.0000 104.0000 105.0000 106.0000 107.0000 108.0000 109.0000 110.0000 111.0000 112.000 113.0000 114.0000 115.0000 116.0000 117.0000 118.0000 119.0000 120.0000 121.0000 122.0000 123.0000 124.0000 125.0000 126.0000 127.0000 128.0000 129.0000 130.0000
## Effort 461.7988 466.4634 471.1280 475.7927 480.4573 485.1220 489.7866 494.4512 499.1159 503.7805 508.4451 513.1098 517.7744 522.439 527.1037 531.7683 536.4329 541.0976 545.7622 550.4268 555.0915 559.7561 564.4207 569.0854 573.7500 578.4146 583.0793 587.7439 592.4085 597.0732 601.7378 606.4024
## Richness 95.9083 96.0591 96.1962 96.4262 96.6301 96.7518 96.8709 97.0113 97.1618 97.3120 97.4559 97.5820 97.6625 97.774 97.8875 98.0212 98.1916 98.3062 98.4148 98.5660 98.6982 98.8264 98.9860 99.0901 99.2508 99.4181 99.5826 99.6636 99.7552 99.8501 99.9673 100.0626
## sd 2.3710 2.3421 2.3244 2.2608 2.1510 2.1341 2.1339 2.1036 2.0830 2.1229 2.1097 2.0828 2.0669 2.041 2.0307 2.0210 2.0420 2.0676 2.0742 2.0749 2.0260 1.9965 2.0042 1.9885 1.9345 1.9464 1.8630 1.8210 1.7939 1.7916 1.7548 1.7520
##
## Sites 131.0000 132.0000 133.0000 134.0000 135.0000 136.0000 137.0000 138.0000 139.0000 140.0000 141.0000 142.0000 143.0000 144.0000 145.0000 146.0000 147.0000 148.0000 149.0000 150.0000 151.0000 152.0000 153.0000 154.0000 155.0000 156.0000 157.0000 158.0000 159.0000 160.0000 161.0000 162.0000
## Effort 611.0671 615.7317 620.3963 625.0610 629.7256 634.3902 639.0549 643.7195 648.3841 653.0488 657.7134 662.3780 667.0427 671.7073 676.3720 681.0366 685.7012 690.3659 695.0305 699.6951 704.3598 709.0244 713.6890 718.3537 723.0183 727.6829 732.3476 737.0122 741.6768 746.3415 751.0061 755.6707
## Richness 100.1148 100.2055 100.3376 100.4429 100.5246 100.6333 100.7329 100.7909 100.8665 100.9216 100.9800 101.0768 101.1577 101.2530 101.3564 101.4664 101.5507 101.6486 101.7269 101.7910 101.8817 101.9670 102.0370 102.1180 102.2777 102.3719 102.4443 102.5209 102.5763 102.6302 102.6951 102.7833
## sd 1.7401 1.7160 1.6178 1.6107 1.6278 1.6061 1.6110 1.6030 1.5834 1.5584 1.5354 1.5105 1.4802 1.4457 1.3989 1.3571 1.2741 1.1569 1.1216 1.1036 1.0892 1.0635 1.0139 0.9504 0.8409 0.7648 0.7066 0.6452 0.6152 0.5658 0.5135 0.4589
##
## Sites 163.0000 164
## Effort 760.3354 765
## Richness 102.8792 103
## sd 0.3024 0vegecomm_richness## vibifq vibifqnotrees carex cyper dicot shade shrub hydro svp ap fqi bryo hydro_rel hydro_rel_notrees sens_rel sens_rel_notree tol_rel tol_rel_notree inv small_tree subcanopy canopy S surveycount
## 1002 47.2642730 24.3506371 9 9 30 25 1 9 1 0.05882353 19.470626 0.005272844 0.015818531 0.031686275 0.039624636 0.000313725 0.419138856 0.766065359 0.000000000 0.139364303 0.03419060 0.19225398 22 3
## 1003 76.1688941 52.6722718 4 4 22 22 1 2 1 0.13333333 22.980970 0.009033832 0.018067663 0.066770532 0.488459280 0.051747162 0.110965566 0.042844425 0.000000000 0.048458150 0.03497367 0.20545653 22 6
## 1004 37.8842226 35.1032304 2 2 21 23 1 9 4 0.19047619 15.086592 0.039605806 0.426344850 0.558710387 0.013978520 0.009159187 0.198727955 0.260426207 0.000000000 0.108695652 0.46386133 0.21743769 28 7
## 1005 82.9198687 65.8474332 11 11 47 52 2 18 3 0.08000000 30.512293 0.007461016 0.060135790 0.128541811 0.168544356 0.088512041 0.133204009 0.196852905 0.000000000 0.114361702 0.10777783 0.17574889 30 3
## 1006 72.0717642 64.7405493 3 3 21 23 1 6 2 0.06666667 26.582281 0.008306874 0.110075316 0.254709325 0.259857491 0.046004015 0.123864727 0.068771091 0.000000000 0.132352941 0.22866471 0.17078115 26 7
## 1008 56.3031565 25.1357383 4 4 11 13 0 3 0 0.25000000 16.842857 0.010954302 0.054771512 0.258754528 0.325945265 0.051750906 0.060577292 0.286182508 0.000000000 0.130841121 0.04519758 0.29052172 24 3
## 1009 45.4587252 35.2832254 5 5 36 40 3 10 0 0.05882353 23.372545 0.005486333 0.109214601 0.114944889 0.010972666 0.003849460 0.611835849 0.643860752 0.000000000 0.119873817 0.12339044 0.13820327 25 3
## 1011 13.3041598 14.0109307 6 6 21 14 0 13 0 0.06976744 9.373602 0.005230308 0.020921232 0.023174523 0.005230308 0.005793631 0.755936399 0.831443966 0.117681928 0.162162162 0.00000000 0.24848942 21 3
## 1012 77.4977507 60.9815207 7 7 28 30 1 6 0 0.04000000 25.285714 0.008523696 0.353818616 0.600318164 0.200562564 0.051050906 0.052505967 0.088796761 0.000000000 0.109965636 0.21892384 0.19566287 21 3
## 1013 48.7604648 40.4112843 7 7 15 22 1 3 0 0.04761905 20.329320 0.031235358 0.078400750 0.100126321 0.078713103 0.040688784 0.258784945 0.329698823 0.000000000 0.219512195 0.00000000 0.27540253 14 3
## 1014 24.8850825 21.4048119 4 4 24 23 1 13 1 0.07142857 12.964988 0.002847867 0.076892408 0.157857811 0.000057000 0.000000000 0.440023922 0.488248363 0.000000000 0.071428571 0.00000000 0.35521137 28 6
## 1015 29.1872412 14.1206164 2 2 13 9 0 2 0 0.10000000 12.159200 0.035009919 0.035009919 0.079702444 0.257206208 0.000000000 0.319990664 0.486184910 0.000000000 0.285714286 0.00000000 0.17139260 26 8
## 1016 6.1773014 6.2905703 2 2 24 12 0 13 2 0.21212121 9.043629 0.000000000 0.015707654 0.016215624 0.006703124 0.006919896 0.902659465 0.931850553 0.001709297 0.071428571 0.00000000 0.27012951 44 3
## 1017 57.8356213 34.9038026 1 2 33 27 1 10 0 0.13793103 22.115051 0.096606419 0.111209715 0.225168877 0.252884151 0.136465986 0.219173004 0.239065663 0.000000000 0.157575758 0.03637554 0.20339806 16 3
## 1018 58.2695721 59.3467604 2 2 38 39 1 12 0 0.10256410 24.999623 0.006024546 0.217838023 0.288496100 0.035789380 0.035548534 0.372955159 0.410980547 0.000000000 0.400000000 0.00000000 0.12057408 33 7
## 1019 50.7262363 18.4208465 4 4 18 16 0 6 0 0.13333333 17.709696 0.009402768 0.018993590 0.074653526 0.018993590 0.036957191 0.181708483 0.455497382 0.000000000 0.290322581 0.03699858 0.15836353 32 7
## 1021 49.5166347 40.4169203 5 5 30 31 1 11 4 0.11764706 18.166431 0.004370948 0.409412108 0.575619149 0.006687550 0.000000000 0.181787718 0.255587193 0.000000000 0.090909091 0.46426353 0.18603507 20 3
## 1022 34.0289308 14.9683193 2 2 17 13 0 1 0 0.08333333 14.058212 0.003540110 0.000000000 0.000000000 0.129618918 0.098373506 0.204447118 0.456594160 0.000000000 0.127819549 0.12032384 0.24034522 35 6
## 1024 54.0397256 25.9199503 2 2 26 20 1 5 0 0.00000000 20.080463 0.003951346 0.007902692 0.022591637 0.248934783 0.011295819 0.216513991 0.404333829 0.000000000 0.081081081 0.00000000 0.26441864 25 3
## 1025 64.6754065 45.8101430 4 4 29 33 1 5 0 0.12500000 23.574758 0.009263548 0.034028099 0.086935942 0.139323761 0.071946986 0.152169214 0.365099400 0.000000000 0.157894737 0.10546062 0.18502495 27 7
## 1026 41.4042709 26.7446869 0 0 20 13 1 6 1 0.22222222 14.971083 0.008514512 0.064903074 0.206543967 0.064260469 0.000000000 0.207079362 0.298057260 0.000000000 0.195555556 0.00000000 0.18063116 27 7
## 1028 50.7125368 32.5246403 5 5 23 24 1 10 3 0.04545454 20.396691 0.006966861 0.094439671 0.314004067 0.092984371 0.007875840 0.479521299 0.428280956 0.000000000 0.161616162 0.09295112 0.19941292 18 3
## 1029 0.8577480 0.8595756 0 0 10 8 0 6 0 0.10000000 6.548462 0.000000000 0.014648816 0.014721654 0.000000000 0.000000000 0.975455957 0.975431413 0.941016686 0.190476190 0.00000000 0.13077558 26 3
## 1030 6.6366366 8.4582943 0 0 12 9 0 5 0 0.33333333 8.249579 0.040040040 0.040840841 0.095595127 0.000000000 0.000000000 0.422422422 0.613870665 0.000000000 0.342105263 0.00000000 0.19200706 25 6
## 1031 86.1822440 76.9382667 5 5 36 41 1 3 3 0.02631579 29.552883 0.009376319 0.016455439 0.030452889 0.291728524 0.515934988 0.070900595 0.113858523 0.000000000 0.139130435 0.03919509 0.14176939 21 3
## 1033 90.5469103 87.9129472 1 1 41 40 2 5 3 0.09375000 31.895261 0.009911785 0.158786798 0.252515434 0.294908646 0.421069224 0.020616513 0.017023512 0.000000000 0.031578947 0.16319081 0.26290521 25 7
## 1034 69.7163721 64.8746123 4 4 55 59 1 23 2 0.07692308 29.634900 0.004247006 0.202030069 0.231641995 0.019238937 0.017189326 0.327444152 0.341351772 0.004247006 0.179104478 0.35823996 0.25887917 28 3
## 1035 69.3713786 45.5529917 9 10 23 30 1 14 2 0.02857143 22.248595 0.006892511 0.034600407 0.085434690 0.427542475 0.042717345 0.131922666 0.257665579 0.000000000 0.230303030 0.11035408 0.19782846 21 3
## 1038 66.9286154 28.3226130 0 0 10 11 1 2 2 0.00000000 16.977494 0.028012512 0.056025025 0.214477212 0.653625286 0.107238606 0.022830197 0.085254692 0.000000000 0.275000000 0.00000000 0.28353301 21 3
## 1039 66.3433378 56.3244555 6 6 55 52 1 27 4 0.08064516 26.131183 0.003231696 0.254323741 0.346897844 0.145426342 0.013224112 0.287987246 0.388406862 0.000000000 0.068181818 0.11838253 0.16168238 22 3
## 1040 34.3195861 25.5902373 2 2 18 14 1 12 1 0.15789474 13.926212 0.011650259 0.648531096 0.713675214 0.020679210 0.000000000 0.255257179 0.280897436 0.000000000 0.324324324 0.00000000 0.22559473 22 3
## 1041 72.4508917 52.4095266 4 4 40 41 1 10 7 0.07500000 26.464671 0.051294008 0.037491257 0.090835483 0.137840989 0.109740712 0.145962540 0.240439113 0.000000000 0.041666667 0.10538593 0.20442155 28 7
## 1043 41.4032009 37.9442588 2 2 18 21 0 14 0 0.13636364 14.234383 0.005054079 0.126402507 0.130355468 0.494036187 0.501667883 0.306024462 0.315594704 0.222379460 0.230769231 0.00000000 0.28416918 30 3
## 1044 19.8247978 11.6533864 0 0 9 5 0 3 1 0.33333333 9.250000 0.026954178 0.026954178 0.099601594 0.000539084 0.000000000 0.567115903 0.601593625 0.000000000 0.276190476 0.00000000 0.20367641 19 3
## 1045 71.7057443 60.6749723 10 10 43 53 1 10 5 0.05660377 27.236578 0.013985199 0.028389954 0.042681431 0.117440709 0.034534114 0.130831537 0.186179346 0.000000000 0.224000000 0.11129807 0.28997505 29 3
## 1046 16.7298468 16.7308176 0 0 1 0 0 5 0 0.20000000 6.123724 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.735939192 0.735981896 0.000000000 0.000000000 0.00000000 0.00000000 32 6
## 1047 91.0710699 71.2182695 5 5 45 48 1 10 3 0.11363636 32.173359 0.015856740 0.025709061 0.078284942 0.381978287 0.126118586 0.038669302 0.098435589 0.000000000 0.134883721 0.12313322 0.21540168 18 3
## 1048 61.6108655 37.9946548 4 4 20 24 1 5 3 0.11764706 19.750829 0.008026326 0.052171121 0.120664583 0.544612997 0.238420841 0.194116703 0.395563256 0.000000000 0.061818182 0.13551587 0.20659002 28 3
## 1050 68.2351091 48.2027445 5 5 25 25 1 8 4 0.10526316 23.859095 0.017453022 0.079062191 0.195877775 0.127494328 0.087345056 0.095060795 0.192274431 0.000000000 0.150943396 0.03500834 0.17121292 21 7
## 1051 60.0000000 25.0000000 0 0 5 6 1 2 0 0.00000000 14.000000 0.010141580 0.030123506 0.298804781 0.692840647 0.000000000 0.000000000 0.000000000 0.000000000 0.107954545 0.00000000 0.21661381 27 7
## 1052 21.8469209 19.8484848 0 0 7 5 0 1 0 0.50000000 6.717514 0.000000000 0.005478852 0.009090909 0.000000000 0.000000000 0.139710717 0.231818182 0.000000000 0.135338346 0.00000000 0.40077160 23 6
## 1053 79.1463360 75.8082659 2 2 54 54 1 20 1 0.08620690 31.544663 0.006684045 0.171735401 0.295427542 0.215861239 0.083649534 0.073234855 0.125982139 0.000000000 0.280000000 0.11111865 0.11976083 33 7
## 1054 49.3529079 42.0231102 0 0 25 22 1 8 0 0.10000000 18.007351 0.009740971 0.280661737 0.345487430 0.017046700 0.002997722 0.023183512 0.025540589 0.000000000 0.197674419 0.07054649 0.28838175 32 7
## 1055 5.7153765 4.5936048 1 1 6 6 0 2 0 0.12500000 4.472136 0.002974951 0.000000000 0.000000000 0.000000000 0.000000000 0.845005057 0.977156566 0.000000000 0.046321526 0.00000000 0.24992981 27 6
## 1057 47.0121070 39.0228034 1 1 18 25 1 5 4 0.10526316 18.762424 0.091991619 0.209843103 0.251269812 0.061327746 0.055076189 0.296826289 0.355425005 0.000000000 0.380952381 0.00000000 0.15870464 22 3
## 1058 32.2451969 26.7484609 1 1 20 16 1 9 0 0.20000000 13.139781 0.006768266 0.151157937 0.412142835 0.020304798 0.000000000 0.091913050 0.250238366 0.000000000 0.133333333 0.27578122 0.26959703 15 3
## 1059 7.4596278 7.7087769 4 4 16 8 0 6 0 0.16666667 5.208554 0.000000000 0.004988327 0.005544221 0.000000000 0.000000000 0.819921384 0.911292469 0.000000000 0.000000000 0.00000000 0.12568772 26 3
## 1060 86.0217979 75.2021735 2 2 27 33 1 3 5 0.00000000 29.340949 0.008283327 0.193277620 0.317710067 0.128962191 0.189294684 0.011375769 0.018699507 0.000000000 0.074380165 0.31247163 0.26409609 22 3
## 1061 64.4979646 44.4512500 7 8 24 29 1 8 1 0.03571429 23.590713 0.021376149 0.078379213 0.160279761 0.214331520 0.088008160 0.185806237 0.234251299 0.000000000 0.145038168 0.05658729 0.18345411 19 7
## 1063 0.6882786 0.7016067 0 0 7 3 0 1 0 0.23076923 2.910427 0.003754247 0.000000000 0.000000000 0.000000000 0.000000000 0.969665684 0.988442625 0.000000000 0.000000000 0.00000000 0.00000000 38 7
## 1064 68.9815761 48.7614881 3 3 14 15 1 2 2 0.09090909 20.199962 0.008518371 0.014197285 0.129735340 0.259810324 0.233523612 0.000851837 0.007784120 0.000000000 0.082901554 0.03687149 0.20108903 34 8
## 1065 36.4339155 26.4460066 6 7 41 25 2 30 3 0.04000000 19.398969 0.005679667 0.166641424 0.194026232 0.005679667 0.002204343 0.745708696 0.808736544 0.000000000 0.056547619 0.00000000 0.17096747 25 3
## 1066 71.6540458 68.8886354 4 4 29 31 1 2 3 0.00000000 27.955103 0.034103698 0.038650858 0.098776909 0.097229643 0.178757154 0.051360169 0.104935940 0.000000000 0.140625000 0.11569459 0.25493198 25 3
## 1067 49.0975392 28.8159047 2 2 19 18 1 5 0 0.10000000 17.826752 0.011034686 0.116084894 0.402602373 0.176996359 0.000000000 0.227424872 0.213930348 0.000000000 0.054347826 0.14472548 0.16772444 22 7
## 1068 86.8148642 77.5296486 6 6 51 53 1 17 5 0.14583333 30.137684 0.049954209 0.134412503 0.163033599 0.385931943 0.394362133 0.096958740 0.108948743 0.000000000 0.241379310 0.08966525 0.20454119 23 8
## 1069 49.5071980 51.5128576 5 5 41 32 1 21 1 0.11904762 22.584706 0.002664960 0.200584444 0.310688341 0.072692731 0.106464500 0.204456553 0.298294729 0.000000000 0.152777778 0.23370100 0.12741797 40 3
## 1070 69.8973526 55.2878470 2 2 20 20 2 10 5 0.00000000 23.334524 0.068522032 0.311195446 0.557770430 0.248875536 0.085025980 0.132932743 0.030420406 0.000000000 0.132911392 0.30080543 0.17851444 14 3
## 1071 44.0556446 31.0731281 3 3 28 24 1 11 1 0.09523809 18.577086 0.017223067 0.076585240 0.157366993 0.202371042 0.026542409 0.185664667 0.381148991 0.000000000 0.148148148 0.03828634 0.09913688 30 7
## 1072 22.7373661 21.3626378 5 6 17 17 1 8 0 0.18421053 12.450604 0.011351382 0.281230490 0.383464465 0.000681083 0.000619115 0.394668634 0.390790662 0.000000000 0.340425532 0.00000000 0.39070309 31 7
## 1073 32.0698544 7.3519025 6 6 35 16 1 7 0 0.09375000 15.748536 0.002367527 0.004782405 0.008101442 0.101803661 0.000000000 0.534212738 0.884908725 0.000000000 0.030470914 0.00000000 0.19053064 35 7
## 1074 30.4423741 13.7202290 0 0 26 15 0 1 0 0.20000000 14.411534 0.004714313 0.000000000 0.000000000 0.009428625 0.000000000 0.410522346 0.673015873 0.000000000 0.162500000 0.00000000 0.42606512 46 7
## 1075 30.9225792 19.1232646 4 5 26 20 0 12 1 0.03225807 14.142857 0.024070607 0.064348756 0.151197134 0.024070607 0.056557532 0.546563252 0.624018098 0.000000000 0.242424242 0.00000000 0.30032618 33 6
## 1076 35.1904406 29.1257065 5 6 34 22 1 22 2 0.07142857 17.008401 0.009766367 0.377361562 0.406555154 0.015878485 0.003507305 0.630220937 0.652671544 0.000000000 0.103734440 0.30960972 0.20621328 25 3
## 1077 49.8635649 24.4490187 8 8 27 21 1 11 2 0.12000000 17.649664 0.004873769 0.159640316 0.441313183 0.063359002 0.013473155 0.289136368 0.704982822 0.000000000 0.092807425 0.17718386 0.20787458 17 3
## 1078 26.2286053 22.1832371 4 4 26 24 1 15 0 0.21875000 14.010778 0.002774310 0.107273316 0.116382900 0.124843945 0.135445616 0.561011698 0.608652467 0.123456790 0.210526316 0.08238879 0.08245940 29 3
## 1080 22.1493130 14.1714373 5 6 28 20 0 14 1 0.03333333 15.100260 0.003637488 0.032737392 0.058499805 0.006529291 0.000000000 0.418602113 0.631017897 0.000000000 0.125954198 0.03460696 0.16864491 32 3
## 1081 67.6441016 33.1461639 4 4 23 21 1 8 0 0.33333333 19.776610 0.006473051 0.120830277 0.520156047 0.514866439 0.028422813 0.181633798 0.558424670 0.000000000 0.067524116 0.06885338 0.47382135 26 7
## 1082 52.2306642 42.5636909 5 5 43 40 1 13 4 0.04444444 24.002841 0.018320051 0.061097371 0.107320998 0.022961131 0.016090105 0.440464922 0.467882363 0.000000000 0.099630996 0.00000000 0.27425889 33 7
## 1083 32.2417471 25.9868553 1 2 13 13 0 9 0 0.07142857 12.701706 0.008325124 0.308029582 0.665228335 0.004245813 0.000000000 0.142359618 0.217188062 0.000000000 0.161290323 0.00000000 0.19785790 29 7
## 1086 59.1996968 52.8073944 1 1 32 26 1 18 0 0.03333333 21.144229 0.004861291 0.636018927 0.708483755 0.048612912 0.043321300 0.112619912 0.114620939 0.009722582 0.196721311 0.47819559 0.17963810 34 6
## 1087 69.9907734 38.0826786 2 2 13 13 0 2 1 0.10000000 18.766739 0.010051261 0.000201025 0.001755926 0.563071665 0.089552239 0.007370925 0.064383963 0.000000000 0.132530120 0.00000000 0.32294425 18 3
## 1088 22.6745597 23.3575929 0 0 11 3 0 2 0 0.07692308 10.253048 0.003005651 0.000000000 0.000000000 0.026048972 0.026856175 0.439947101 0.453580135 0.000000000 0.000000000 0.00000000 0.00000000 25 3
## 1089 92.2211128 61.0635030 4 4 44 41 1 9 6 0.19444444 30.250000 0.087474376 0.029423199 0.088148874 0.603467166 0.183286126 0.069493532 0.096699050 0.000000000 0.095238095 0.00000000 0.28637885 26 7
## 1090 55.8924083 31.2538117 3 3 34 29 1 9 0 0.17391304 20.164033 0.017208150 0.027808370 0.095081196 0.323788546 0.047540598 0.175201909 0.551973013 0.000000000 0.208695652 0.22057155 0.17443767 37 7
## 1091 55.3715779 33.3228587 6 6 13 15 1 6 1 0.00000000 14.818657 0.015372134 0.208411973 0.718504343 0.389529869 0.017665244 0.031359153 0.107757986 0.000000000 0.130434783 0.14917110 0.20109531 19 3
## 1092 63.1114138 19.1674678 3 3 15 12 0 1 0 0.00000000 18.557783 0.019445457 0.007058242 0.042725913 0.342395332 0.064516129 0.077922995 0.471694082 0.000000000 0.156950673 0.03491555 0.18949034 19 3
## 1093 56.1844921 26.1134860 3 3 19 13 0 4 1 0.09090909 18.278154 0.009517766 0.019225888 0.179396092 0.095748731 0.092362345 0.191878173 0.454706927 0.000000000 0.272727273 0.00000000 0.15406053 23 3
## 1094 34.6268174 21.0755361 1 1 33 24 1 9 1 0.10000000 17.259516 0.010561521 0.031790178 0.081041059 0.018641084 0.000000000 0.393416652 0.437514023 0.000000000 0.086419753 0.03497019 0.17006012 31 3
## 1097 75.3499787 45.2044752 2 2 27 21 0 4 1 0.10526316 24.833333 0.005980980 0.017942941 0.131262306 0.495015850 0.055348939 0.041986483 0.307153796 0.000000000 0.078804348 0.13952624 0.16486691 33 7
## 1098 65.0847028 62.3960785 4 4 58 50 1 21 2 0.06349206 26.052335 0.002918827 0.171354628 0.190698933 0.003093957 0.003378269 0.240355708 0.228444589 0.002918827 0.162790698 0.03799072 0.13694672 41 7
## 1099 57.6102830 20.0327567 3 3 15 11 0 1 0 0.12500000 16.341190 0.010840892 0.000000000 0.000000000 0.553319120 0.057164207 0.156614751 0.417336073 0.000000000 0.201550388 0.03598872 0.22731842 24 3
## 1100 28.3592783 7.1167580 1 1 12 7 0 5 1 0.00000000 11.783766 0.014668492 0.029336984 0.074738416 0.132309799 0.000747384 0.514668492 0.861235675 0.000000000 0.236714976 0.00000000 0.19891336 20 3
## 1101 66.6230627 53.7345943 5 5 43 39 0 0 2 0.04761905 27.287642 0.002503703 0.000000000 0.000000000 0.047057105 0.025676038 0.309232406 0.358469679 0.000000000 0.013392857 0.17047783 0.20196864 24 7
## 1102 43.5536105 26.8481837 5 5 21 27 0 0 0 0.05555556 16.302936 0.011794926 0.000000000 0.000000000 0.009304886 0.024861685 0.179544978 0.479305273 0.007863284 0.220512821 0.00000000 0.14265512 33 6
## 1103 51.0794746 23.7739794 8 9 22 23 0 0 2 0.06250000 18.243355 0.007918910 0.000000000 0.000000000 0.079347482 0.042271382 0.219934537 0.587008595 0.000000000 0.131034483 0.03510817 0.21502618 12 1
## 1104 17.4910984 14.0366974 6 6 31 22 0 20 1 0.08571429 15.105187 0.002750313 0.027503128 0.029980063 0.005500626 0.002998006 0.841664489 0.863500772 0.000000000 0.060498221 0.00000000 0.27189190 30 3
## 1105 64.6964124 54.1651226 8 8 27 33 0 0 1 0.07692308 24.651184 0.015704340 0.000000000 0.000000000 0.171988693 0.170685279 0.206040936 0.177203507 0.000000000 0.222222222 0.06094992 0.18350302 29 7
## 1107 7.5175044 6.5512624 0 0 12 9 0 0 0 0.94736842 4.585303 0.000000000 0.000000000 0.000000000 0.018754689 0.000000000 0.772693173 0.788064269 0.018754689 0.000000000 0.00000000 0.00000000 25 3
## 1108 39.3594051 36.4378180 2 2 34 22 2 39 1 0.13953488 20.579116 0.000000000 0.101119398 0.103930836 0.072136131 0.070169864 0.775862735 0.797434166 0.000000000 0.000000000 0.16789173 0.30721036 32 3
## 1110 50.3361915 38.9448141 3 3 16 17 0 0 0 0.07692308 16.614251 0.012599483 0.000000000 0.000000000 0.021083136 0.000000000 0.056277693 0.138987289 0.000000000 0.116666667 0.42389991 0.21260174 20 3
## 1111 40.3099379 36.7915692 2 2 21 23 0 0 1 0.05555556 18.584678 0.019406171 0.000000000 0.000000000 0.077624685 0.079978672 0.349311081 0.399893362 0.000000000 0.233766234 0.00000000 0.21955821 19 3
## 1112 49.5745936 26.5228283 6 6 30 28 1 15 2 0.05405405 19.564659 0.032777482 0.098496333 0.198590153 0.180276150 0.049564930 0.286283990 0.461559643 0.000000000 0.123076923 0.19252234 0.21024158 26 3
## 1113 67.3125523 44.0752952 1 1 19 17 1 1 0 0.11111111 22.299904 0.006195275 0.004171485 0.038055765 0.260325458 0.169555388 0.241739633 0.170685757 0.000000000 0.052380952 0.04775089 0.24970445 27 7
## 1114 38.5917371 14.6986197 2 2 16 14 0 0 1 0.16666667 15.493175 0.012769354 0.000000000 0.000000000 0.221744736 0.000000000 0.241696851 0.469809069 0.000000000 0.203883495 0.00000000 0.17589600 21 7
## 1115 37.0430682 16.6196300 1 1 11 10 0 0 0 0.10000000 13.764944 0.033090668 0.000000000 0.000000000 0.083057578 0.047977422 0.036399735 0.101599247 0.000000000 0.068493151 0.00000000 0.22685374 17 3
## 1116 22.0782247 22.8754990 5 5 27 21 0 0 2 0.18181818 13.567731 0.003262536 0.000000000 0.000000000 0.009787609 0.010380623 0.613356824 0.650519031 0.000000000 0.000000000 0.00000000 0.00000000 33 6
## 1117 80.3722481 54.2961888 3 3 24 30 0 0 2 0.14814815 24.363107 0.008905116 0.000000000 0.000000000 0.445612004 0.105351171 0.023420455 0.109113712 0.000000000 0.047732697 0.03492895 0.21164164 24 7
## 1119 5.0552703 4.7242690 2 2 15 4 1 3 0 0.07692308 8.966120 0.000000000 0.000000000 0.000000000 0.011494693 0.000000000 0.889516840 0.968385919 0.000000000 0.018867925 0.00000000 0.15788103 27 6
## 1122 49.1094343 36.7583160 9 10 22 25 0 0 0 0.02631579 20.178224 0.012632430 0.000000000 0.000000000 0.071667986 0.051813770 0.322564889 0.441015780 0.012632430 0.383561644 0.00000000 0.18456142 30 6
## 1123 3.9409783 3.9054676 4 5 16 10 0 0 1 0.12500000 8.198916 0.000000000 0.000000000 0.000000000 0.001956954 0.001963698 0.974999918 0.978360051 0.013698675 0.000000000 0.00000000 0.00000000 36 3
## 1124 26.7138853 6.3973468 3 4 26 19 0 10 2 0.00000000 14.320780 0.005453703 0.016579256 0.031289133 0.223710881 0.000205850 0.485470427 0.878394373 0.000000000 0.079787234 0.03600028 0.11806617 28 3
## 1125 44.0553690 32.2007370 7 7 30 32 0 0 3 0.09375000 20.445485 0.012575242 0.000000000 0.000000000 0.041917473 0.000000000 0.448282222 0.557006820 0.000000000 0.136842105 0.03690905 0.14177715 26 3
## 1126 26.0212404 16.7881897 4 4 23 17 0 0 1 0.10000000 14.142857 0.017190992 0.000000000 0.000000000 0.086642599 0.045945325 0.530256146 0.593728463 0.000000000 0.315217391 0.00000000 0.13719234 32 6
## 1127 68.5741717 41.9298087 1 1 8 11 0 0 2 0.00000000 18.432707 0.016758279 0.000000000 0.000000000 0.643718997 0.139013008 0.015384100 0.045576408 0.000000000 0.118644068 0.00000000 0.43130675 37 8
## 1128 64.7408687 52.9120284 2 2 23 25 1 3 2 0.09090909 23.509055 0.019483682 0.055528495 0.085675635 0.272771554 0.210431384 0.166975158 0.257628138 0.000000000 0.063492063 0.13589695 0.26069283 19 6
## 1129 74.6561865 66.0798092 0 0 15 18 1 1 1 0.09090909 23.337820 0.011290717 0.000225814 0.000634350 0.367625750 0.667336258 0.045388683 0.064069356 0.000000000 0.056074766 0.08701678 0.20152191 27 7
## 1131 45.5803365 25.6452461 1 1 12 10 0 0 3 0.09090909 15.652476 0.026589852 0.000000000 0.000000000 0.020607135 0.000000000 0.080833149 0.188819876 0.000000000 0.081395349 0.00000000 0.29849496 28 3
## 1132 61.5461263 48.7898138 0 0 29 32 0 0 3 0.03571429 23.193102 0.008254912 0.000000000 0.000000000 0.181608057 0.068073519 0.227587915 0.301452235 0.000000000 0.205882353 0.44181765 0.28576944 23 6
## 1133 11.0775751 11.0812355 2 3 17 5 0 0 0 0.07407407 6.761234 0.000000000 0.000000000 0.000000000 0.001251951 0.001255961 0.778054134 0.780520889 0.000000000 0.000000000 0.00000000 0.00000000 31 3
## 1134 6.3143417 4.1146735 3 3 17 12 0 0 0 0.10526316 9.328719 0.010534633 0.000000000 0.000000000 0.010534633 0.000000000 0.758581336 0.867277354 0.000000000 0.224719101 0.00000000 0.20337584 31 7
## 1135 62.4428095 56.3290470 7 7 37 34 0 0 0 0.05263158 22.648834 0.005368983 0.000000000 0.000000000 0.025270015 0.012427073 0.162375954 0.149145382 0.000000000 0.108391608 0.24237264 0.18391901 30 7
## 1136 48.4135472 26.0122411 1 1 10 4 0 0 0 0.00000000 15.588457 0.018674718 0.000000000 0.000000000 0.467241433 0.141242938 0.299168975 0.002824859 0.000000000 0.312169312 0.03520657 0.17819956 36 7
## 1138 14.4219204 11.3714328 3 3 35 20 0 0 0 0.07317073 11.324417 0.002655854 0.000000000 0.000000000 0.002655854 0.002834226 0.704624728 0.751948530 0.000000000 0.003267974 0.00000000 0.20428325 38 7
## 1140 48.2773792 35.8349200 3 3 11 15 1 1 2 0.07692308 17.800000 0.024457851 0.024457851 0.041006014 0.305723137 0.164024057 0.082830589 0.138873701 0.000000000 0.207207207 0.00000000 0.30972068 9 2
## 1141 19.2394574 19.7753194 7 9 12 15 0 0 2 0.07142857 11.833333 0.005878780 0.000000000 0.000000000 0.017636339 0.018281164 0.832337207 0.862769393 0.005878780 0.416243655 0.07767635 0.06882081 35 4
## 1143 56.2833929 44.5067290 3 3 21 26 0 0 2 0.04761905 22.062113 0.013690191 0.000000000 0.000000000 0.082414950 0.059209881 0.113697036 0.162827174 0.000000000 0.157894737 0.03578353 0.22902443 24 3
## 1144 64.4036676 50.7355473 7 7 29 36 1 12 3 0.02631579 24.270696 0.009465215 0.132513015 0.243013366 0.104306673 0.052074293 0.223284430 0.270612741 0.000000000 0.046920821 0.12178217 0.15508403 18 2
## 1145 74.0871563 55.5889607 2 2 29 29 1 9 3 0.07692308 24.675454 0.027309516 0.143052211 0.376281591 0.227438616 0.082348331 0.062563619 0.164566055 0.000000000 0.013245033 0.08775583 0.30325680 27 7
## 1146 48.7429834 33.8597926 7 7 23 25 0 0 2 0.00000000 20.833333 0.042099354 0.000000000 0.000000000 0.042099354 0.000000000 0.151557676 0.179257362 0.000000000 0.285714286 0.00000000 0.17181709 20 3
## 1147 55.8073915 35.4732695 3 3 12 14 0 0 1 0.11111111 18.371173 0.022140221 0.000000000 0.000000000 0.067306273 0.057475582 0.045166052 0.114951165 0.000000000 0.094339623 0.07925455 0.15954906 26 7
## 1148 31.3509544 11.9440054 0 0 19 10 0 3 0 0.00000000 12.247449 0.012371644 0.000000000 0.000000000 0.024743288 0.000000000 0.435399398 0.644840822 0.000000000 0.162162162 0.14737364 0.23031636 20 3
## 1149 31.6438090 14.7058823 2 2 15 12 0 0 0 0.14285714 13.761951 0.018796992 0.000000000 0.000000000 0.125563910 0.000000000 0.206766917 0.294117647 0.000000000 0.171171171 0.00000000 0.23571960 39 7
## 1150 13.7657362 9.2470760 0 0 7 5 0 0 0 0.00000000 6.948792 0.082101806 0.000000000 0.000000000 0.082101806 0.000000000 0.415435140 0.554824561 0.000000000 0.562500000 0.00000000 0.17073555 29 6
## 1151 53.0986317 19.0950517 3 3 12 10 0 0 1 0.00000000 15.600270 0.012645955 0.000000000 0.000000000 0.089027526 0.004713276 0.038443704 0.358208955 0.000000000 0.118421053 0.00000000 0.17254395 20 3
## 1153 64.9326107 51.9177502 2 2 27 29 0 0 3 0.12500000 24.019223 0.020253848 0.000000000 0.000000000 0.175533351 0.149812734 0.108290575 0.200249688 0.000000000 0.103448276 0.00000000 0.42040109 29 7
## 1154 55.5187565 45.7347820 0 0 21 21 0 0 0 0.05263158 16.537377 0.004964504 0.000000000 0.000000000 0.084396565 0.021809062 0.047659236 0.052341748 0.000000000 0.022421525 0.36813310 0.13428338 28 3
## 1155 40.3122142 24.3452125 2 2 12 15 0 0 1 0.07142857 12.894156 0.012861185 0.000000000 0.000000000 0.285175341 0.000000000 0.142244705 0.237440962 0.000000000 0.137254902 0.00000000 0.23696839 24 3
## 1156 1.0813491 1.0265884 1 1 17 14 0 12 0 0.09523809 9.128709 0.002620923 0.015725537 0.015892145 0.000000000 0.000000000 0.952770971 0.962865354 0.002620923 0.478260870 0.00000000 0.24308641 30 3
## 1157 61.5372140 40.1781105 7 7 15 13 0 0 0 0.11764706 19.843135 0.010522439 0.000000000 0.000000000 0.226653338 0.142510121 0.081583978 0.192442645 0.000000000 0.071770335 0.03514090 0.19232775 17 3
## 1159 34.1045005 20.7453541 6 6 23 22 0 0 1 0.06896552 15.776277 0.029364551 0.000000000 0.000000000 0.047218198 0.010911252 0.367879096 0.466861125 0.000000000 0.217777778 0.10655868 0.19307106 16 1
## 1161 44.2238010 36.0447104 2 2 14 14 1 7 3 0.10000000 16.819419 0.027068483 0.135883786 0.262827225 0.198502211 0.104712042 0.231976902 0.291623037 0.000000000 0.240000000 0.00000000 0.19752813 29 7
## 1162 84.7094202 82.8188118 7 8 61 49 0 0 4 0.09090909 28.781504 0.006421576 0.000000000 0.000000000 0.404719859 0.477399404 0.284582865 0.281107633 0.000000000 0.166666667 0.11502740 0.16604261 27 3
## 1164 27.9469239 10.6513497 2 2 23 16 0 7 0 0.00000000 12.809880 0.009602612 0.016260423 0.025736505 0.177936399 0.000000000 0.418401805 0.601134838 0.009602612 0.061538462 0.00000000 0.20613724 24 3
## 1166 46.4687101 42.2716973 2 2 20 22 0 0 2 0.04166667 19.438548 0.021276596 0.000000000 0.000000000 0.063829787 0.051948052 0.212765957 0.259740260 0.021276596 0.269230769 0.00000000 0.26124182 28 3
## 1167 15.6939345 9.4909865 4 5 24 17 0 0 1 0.02564103 13.483997 0.001710888 0.000000000 0.000000000 0.005166883 0.003714768 0.861229844 0.934935843 0.000000000 0.055016181 0.08762413 0.17993393 20 3
## 1168 33.4930442 32.1376096 2 2 25 25 1 12 0 0.02941176 16.130857 0.002667876 0.262172185 0.265014428 0.005335752 0.005393598 0.539617960 0.545468029 0.002667876 0.046511628 0.08721168 0.21476748 25 3
## 1169 44.8892920 35.6223573 5 5 24 27 0 0 2 0.11111111 19.975541 0.016398371 0.000000000 0.000000000 0.049523081 0.045172219 0.311897018 0.339243365 0.000000000 0.277777778 0.00000000 0.23198093 20 3
## 1171 51.1157162 32.4117032 2 2 21 26 0 0 1 0.13636364 18.493242 0.007276167 0.000000000 0.000000000 0.109288027 0.094694479 0.146360098 0.317463240 0.000000000 0.164473684 0.18532007 0.21561005 24 3
## 1173 11.7999514 4.9824524 1 1 32 21 0 0 3 0.20408163 13.076697 0.002379272 0.000000000 0.000000000 0.009517087 0.000000000 0.873704288 0.927745079 0.130859948 0.120879121 0.00000000 0.17192061 27 3
## 1174 45.4563394 41.3060154 6 8 38 31 0 0 2 0.10204082 17.094791 0.004708172 0.000000000 0.000000000 0.009416344 0.005153311 0.355466972 0.353001804 0.000000000 0.063492063 0.55001752 0.21682028 22 3
## 1175 47.1780885 30.1299671 2 2 26 24 0 0 0 0.20689655 19.990242 0.020109934 0.000000000 0.000000000 0.127631050 0.087834870 0.228448854 0.374176548 0.000000000 0.082568807 0.03527536 0.20087728 18 3
## 1176 79.3295904 62.5540727 5 5 22 26 1 3 4 0.05555556 27.666667 0.023475097 0.017723698 0.040644206 0.361105677 0.134852631 0.039281662 0.049706159 0.000000000 0.082758621 0.03568503 0.18198608 25 6
## 1177 10.8313493 4.2607287 1 1 21 16 0 7 1 0.00000000 11.672328 0.002994610 0.013475744 0.015087170 0.020962268 0.000000000 0.659612697 0.735136343 0.000000000 0.036659878 0.00000000 0.24035969 21 3
## 1178 63.5212313 48.0175070 5 5 15 20 0 0 3 0.00000000 23.186922 0.017584994 0.000000000 0.000000000 0.211019930 0.100446429 0.070691676 0.067633929 0.000000000 0.094936709 0.12220624 0.18210034 29 7
## 1179 19.2360846 11.6142806 3 3 24 19 0 0 0 0.10810811 13.152186 0.000000000 0.000000000 0.000000000 0.009867123 0.008879242 0.470727536 0.552959747 0.002466781 0.093525180 0.00000000 0.25779149 26 3
## 1180 1.8083853 2.5876619 4 4 13 13 0 0 0 0.00000000 9.449112 0.003971564 0.000000000 0.000000000 0.003971564 0.005689253 0.644386195 0.923081299 0.003971564 0.062015504 0.00000000 0.25432460 30 6
## 1182 47.4456174 32.5749792 0 0 16 13 0 0 0 0.11111111 17.000000 0.019193858 0.000000000 0.000000000 0.038771593 0.063511831 0.192706334 0.250311333 0.000000000 0.283185841 0.00000000 0.18594614 32 7
## 1183 22.1725624 15.4762436 2 2 16 11 1 4 0 0.00000000 10.253048 0.014928344 0.044785032 0.080303368 0.014928344 0.000000000 0.293665406 0.526031675 0.000000000 0.212765957 0.00000000 0.26602949 23 6
## 1184 36.0002961 20.7327971 6 6 10 7 0 6 0 0.16666667 12.578842 0.028764562 0.014669927 0.024859859 0.273263340 0.000000000 0.179778513 0.255910310 0.000000000 0.105485232 0.00000000 0.23273561 16 3
## 1187 31.4765682 24.5520339 2 2 28 24 0 0 2 0.14285714 17.195706 0.008974915 0.000000000 0.000000000 0.049362033 0.042464562 0.511570161 0.655167533 0.000000000 0.017953321 0.00000000 0.20523612 22 3
## 1188 11.4124737 11.6851910 1 1 22 11 0 9 1 0.06451613 9.170605 0.003078328 0.328416559 0.343207639 0.001067154 0.001115216 0.918732139 0.934373844 0.000000000 0.020408163 0.00000000 0.09499185 27 3
## 1189 78.2998025 69.5002046 12 12 49 59 0 0 7 0.08571429 31.045780 0.006506632 0.000000000 0.000000000 0.176980402 0.089774282 0.233436283 0.265128035 0.000000000 0.188888889 0.15676817 0.13104936 30 3
## 1190 23.0916312 21.6842615 5 5 26 29 0 0 0 0.10256410 16.424201 0.004528097 0.000000000 0.000000000 0.031696678 0.032531408 0.754682807 0.774557340 0.566012105 0.166666667 0.05848294 0.12946399 32 6
## 1191 48.9535176 37.1858220 5 5 24 26 0 0 2 0.04166667 20.739806 0.012896570 0.000000000 0.000000000 0.077379417 0.046468401 0.202390164 0.364622057 0.000000000 0.059090909 0.07114911 0.11167694 28 8
## 1192 72.4031045 50.2791237 4 5 37 32 1 12 5 0.07894737 23.576676 0.014490654 0.294643288 0.532611543 0.403323190 0.034925347 0.103994590 0.187985681 0.000000000 0.014150943 0.22526683 0.32548705 30 6
## 1194 21.3756974 8.5813381 2 2 27 17 0 0 2 0.05555556 14.388363 0.006841505 0.000000000 0.000000000 0.006841505 0.000000000 0.464082098 0.564371257 0.000000000 0.065868263 0.08676337 0.22470780 21 3
## 1195 24.4356468 6.7525974 3 4 30 24 0 0 0 0.10526316 14.293825 0.004343231 0.000000000 0.000000000 0.078178159 0.000000000 0.616912542 0.829813786 0.000000000 0.047272727 0.00000000 0.20777354 27 7
## 1196 33.4535675 23.5771189 6 6 34 29 0 0 0 0.15384615 16.016063 0.006009435 0.000000000 0.000000000 0.186292479 0.054861657 0.497631281 0.656017409 0.000000000 0.144067797 0.00000000 0.22926804 35 8
## 1197 74.0828500 58.7189028 1 1 27 34 0 0 4 0.07142857 26.231569 0.020425764 0.000000000 0.000000000 0.272479688 0.092109303 0.082247742 0.139085048 0.000000000 0.049180328 0.25528860 0.26281393 24 3
## 1198 10.5315333 9.5710090 3 4 19 11 0 0 0 0.00000000 10.224154 0.012729124 0.000000000 0.000000000 0.000000000 0.000000000 0.681262729 0.763373191 0.000000000 0.129921260 0.00000000 0.27265831 20 3
## 1199 72.7472391 68.0281154 8 8 27 32 0 0 2 0.09677419 24.917691 0.042283960 0.000000000 0.000000000 0.260218624 0.279594541 0.156654242 0.179143236 0.000000000 0.323809524 0.03651024 0.13588334 19 3
## 1200 51.7404400 39.5929323 1 1 9 9 0 0 0 0.25000000 15.377079 0.130737943 0.000000000 0.000000000 0.581057525 0.489457831 0.022225450 0.001129518 0.000000000 0.155555556 0.18603820 0.23530338 32 6
## 1257 62.8582229 36.3095238 1 1 14 12 1 1 2 0.09090909 21.059956 0.025000000 0.025000000 0.047619048 0.437500000 0.190476190 0.125000000 0.238095238 0.000000000 0.188888889 0.00000000 0.36842214 32 3
## 3980 1.1376107 1.0261715 0 0 4 3 0 0 0 0.00000000 2.474874 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.951658235 0.960496503 0.000000000 0.086956522 0.00000000 0.21715732 23 6CMP_bbs_pcap_join_pcaponly_pcapsum_date<-CMP_bbs_pcap_join_pcaponly %>%
group_by(simp_plot, reservatio,date,species_4_) %>%
summarise(n=n(),longitude=mean(xcoord), latitude=mean(ycoord),vibifq=mean(vibi_fq),vibifqnotrees=mean(vibi_fq_no_trees),carex=mean(carex_metric_value),cyper=mean(cyperaceae_metric_value),
dicot=mean(dicot_metric_value),shade=mean(shade_metric_value),shrub=mean(shrub_metric_value),
hydro=mean(hydrophyte_metric_value),svp=mean(svp_metric_value),ap=mean(ap_ratio_metric_value),
fqi=mean(fqai_score),bryo=mean(bryophyte_metric_value),hydro_rel=mean(hydrophyte_rel_cov_metric_value),
hydro_rel_notrees=mean(hydrophyte_rel_cov_no_trees),
sens_rel=mean(sensitive_rel_cov_metric_value),
sens_rel_notree=mean(sensitive_rel_cov_no_trees),
tol_rel=mean(tolerant_rel_cov_metric_value),
tol_rel_notree=mean(tolerant_rel_cov_no_trees),
inv=mean(invasive_graminoids_metric_value),
small_tree=mean(small_tree_metric_value),
subcanopy=mean(subcanopy_iv),
canopy=mean(canopy_iv.1)
) %>%
spread(species_4_,n) %>%
replace(is.na(.),0)
birds_bycount<-CMP_bbs_pcap_join_pcaponly_pcapsum_date[,29:ncol(CMP_bbs_pcap_join_pcaponly_pcapsum_date)]
vegecomm_bycount<-CMP_bbs_pcap_join_pcaponly_pcapsum_date[,1:28]
index<-as.numeric(factor(CMP_bbs_pcap_join_pcaponly_pcapsum_date$simp_plot))
index_unique<-unique(index)
rarefylist<-list()
for(i in index_unique){
birds_temp<-birds_bycount[index %in% i,]
rarefylist[[i]]<- specaccum(birds_temp)
}
totalspecaccum<-specaccum(birds_bycount)
plot(totalspecaccum)
plot(rarefylist[[1]])
So the number of surveys clearly has an effect on the species richness of the site. This is somewhat difficult to account for, but we will attempt to do so using mixed effects models. I then use Akaike’s Information Criterion (AIC for short) to rank these models. The lowest number is the best, and any model within 2 AIC of the top model is potentially worth interpreting.
simp_plot<-factor(CMP_bbs_pcap_join_pcaponly_pcapsum_date$simp_plot)
date<-strptime(CMP_bbs_pcap_join_pcaponly_pcapsum_date$date,format="%m/%d/%Y")
month<-month(date)
S_plotdate<-specnumber(birds_bycount)
Rich_vibifq<-lmer(S_plotdate~vibifq+(1|simp_plot),vegecomm_bycount)
summary(Rich_vibifq)## Linear mixed model fit by REML ['lmerMod']
## Formula: S_plotdate ~ vibifq + (1 | simp_plot)
## Data: vegecomm_bycount
##
## REML criterion at convergence: 4350.7
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -2.6215 -0.6926 -0.1170 0.6499 3.4453
##
## Random effects:
## Groups Name Variance Std.Dev.
## simp_plot (Intercept) 5.692 2.386
## Residual 13.710 3.703
## Number of obs: 765, groups: simp_plot, 164
##
## Fixed effects:
## Estimate Std. Error t value
## (Intercept) 14.42728 0.51881 27.809
## vibifq -0.02997 0.01009 -2.972
##
## Correlation of Fixed Effects:
## (Intr)
## vibifq -0.892Rich_vibifqnotrees<-lmer(S_plotdate~vibifqnotrees+(1|simp_plot),vegecomm_bycount)
summary(Rich_vibifqnotrees)## Linear mixed model fit by REML ['lmerMod']
## Formula: S_plotdate ~ vibifqnotrees + (1 | simp_plot)
## Data: vegecomm_bycount
##
## REML criterion at convergence: 4356.9
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -2.5584 -0.6948 -0.1026 0.6376 3.4659
##
## Random effects:
## Groups Name Variance Std.Dev.
## simp_plot (Intercept) 6.088 2.467
## Residual 13.702 3.702
## Number of obs: 765, groups: simp_plot, 164
##
## Fixed effects:
## Estimate Std. Error t value
## (Intercept) 13.61824 0.46208 29.47
## vibifqnotrees -0.01692 0.01192 -1.42
##
## Correlation of Fixed Effects:
## (Intr)
## vibifqnotrs -0.855Rich_carex<-lmer(S_plotdate~carex+(1|simp_plot),vegecomm_bycount)
summary(Rich_carex)## Linear mixed model fit by REML ['lmerMod']
## Formula: S_plotdate ~ carex + (1 | simp_plot)
## Data: vegecomm_bycount
##
## REML criterion at convergence: 4353.8
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -2.5078 -0.6891 -0.0956 0.6328 3.5149
##
## Random effects:
## Groups Name Variance Std.Dev.
## simp_plot (Intercept) 6.196 2.489
## Residual 13.685 3.699
## Number of obs: 765, groups: simp_plot, 164
##
## Fixed effects:
## Estimate Std. Error t value
## (Intercept) 12.74937 0.40130 31.770
## carex 0.09391 0.09746 0.964
##
## Correlation of Fixed Effects:
## (Intr)
## carex -0.800Rich_cyper<-lmer(S_plotdate~cyper+(1|simp_plot),vegecomm_bycount)
summary(Rich_cyper)## Linear mixed model fit by REML ['lmerMod']
## Formula: S_plotdate ~ cyper + (1 | simp_plot)
## Data: vegecomm_bycount
##
## REML criterion at convergence: 4353.7
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -2.5164 -0.6870 -0.0952 0.6341 3.5192
##
## Random effects:
## Groups Name Variance Std.Dev.
## simp_plot (Intercept) 6.186 2.487
## Residual 13.685 3.699
## Number of obs: 765, groups: simp_plot, 164
##
## Fixed effects:
## Estimate Std. Error t value
## (Intercept) 12.71752 0.39948 31.84
## cyper 0.09949 0.09302 1.07
##
## Correlation of Fixed Effects:
## (Intr)
## cyper -0.798Rich_dicot<-lmer(S_plotdate~dicot+(1|simp_plot),vegecomm_bycount)
summary(Rich_dicot)## Linear mixed model fit by REML ['lmerMod']
## Formula: S_plotdate ~ dicot + (1 | simp_plot)
## Data: vegecomm_bycount
##
## REML criterion at convergence: 4351.4
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -2.4704 -0.7081 -0.0930 0.6339 3.5457
##
## Random effects:
## Groups Name Variance Std.Dev.
## simp_plot (Intercept) 5.934 2.436
## Residual 13.667 3.697
## Number of obs: 765, groups: simp_plot, 164
##
## Fixed effects:
## Estimate Std. Error t value
## (Intercept) 11.80814 0.54551 21.646
## dicot 0.05250 0.02068 2.539
##
## Correlation of Fixed Effects:
## (Intr)
## dicot -0.900Rich_shade<-lmer(S_plotdate~shade+(1|simp_plot),vegecomm_bycount)
summary(Rich_shade)## Linear mixed model fit by REML ['lmerMod']
## Formula: S_plotdate ~ shade + (1 | simp_plot)
## Data: vegecomm_bycount
##
## REML criterion at convergence: 4356.2
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -2.4711 -0.6893 -0.0932 0.6487 3.5347
##
## Random effects:
## Groups Name Variance Std.Dev.
## simp_plot (Intercept) 6.176 2.485
## Residual 13.676 3.698
## Number of obs: 765, groups: simp_plot, 164
##
## Fixed effects:
## Estimate Std. Error t value
## (Intercept) 12.48609 0.49794 25.075
## shade 0.02613 0.01990 1.313
##
## Correlation of Fixed Effects:
## (Intr)
## shade -0.876Rich_hydro<-lmer(S_plotdate~hydro+(1|simp_plot),vegecomm_bycount)
summary(Rich_hydro)## Linear mixed model fit by REML ['lmerMod']
## Formula: S_plotdate ~ hydro + (1 | simp_plot)
## Data: vegecomm_bycount
##
## REML criterion at convergence: 4349.4
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -2.4429 -0.6824 -0.0808 0.6376 3.4771
##
## Random effects:
## Groups Name Variance Std.Dev.
## simp_plot (Intercept) 5.831 2.415
## Residual 13.684 3.699
## Number of obs: 765, groups: simp_plot, 164
##
## Fixed effects:
## Estimate Std. Error t value
## (Intercept) 12.52655 0.30454 41.132
## hydro 0.09721 0.03549 2.739
##
## Correlation of Fixed Effects:
## (Intr)
## hydro -0.632Rich_shrub<-lmer(S_plotdate~shrub+(1|simp_plot),vegecomm_bycount)
summary(Rich_shrub)## Linear mixed model fit by REML ['lmerMod']
## Formula: S_plotdate ~ shrub + (1 | simp_plot)
## Data: vegecomm_bycount
##
## REML criterion at convergence: 4351.3
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -2.5380 -0.6891 -0.0877 0.6446 3.4829
##
## Random effects:
## Groups Name Variance Std.Dev.
## simp_plot (Intercept) 6.204 2.491
## Residual 13.691 3.700
## Number of obs: 765, groups: simp_plot, 164
##
## Fixed effects:
## Estimate Std. Error t value
## (Intercept) 13.2033 0.3119 42.34
## shrub -0.3018 0.4136 -0.73
##
## Correlation of Fixed Effects:
## (Intr)
## shrub -0.635Rich_svp<-lmer(S_plotdate~svp+(1|simp_plot),vegecomm_bycount)
summary(Rich_svp)## Linear mixed model fit by REML ['lmerMod']
## Formula: S_plotdate ~ svp + (1 | simp_plot)
## Data: vegecomm_bycount
##
## REML criterion at convergence: 4353.8
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -2.5198 -0.6806 -0.1005 0.6550 3.4904
##
## Random effects:
## Groups Name Variance Std.Dev.
## simp_plot (Intercept) 6.248 2.500
## Residual 13.686 3.699
## Number of obs: 765, groups: simp_plot, 164
##
## Fixed effects:
## Estimate Std. Error t value
## (Intercept) 13.02335 0.31742 41.029
## svp 0.02676 0.15316 0.175
##
## Correlation of Fixed Effects:
## (Intr)
## svp -0.649Rich_ap<-lmer(S_plotdate~ap+(1|simp_plot),vegecomm_bycount)
summary(Rich_ap)## Linear mixed model fit by REML ['lmerMod']
## Formula: S_plotdate ~ ap + (1 | simp_plot)
## Data: vegecomm_bycount
##
## REML criterion at convergence: 4348.3
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -2.5238 -0.6858 -0.0978 0.6528 3.4990
##
## Random effects:
## Groups Name Variance Std.Dev.
## simp_plot (Intercept) 6.234 2.497
## Residual 13.691 3.700
## Number of obs: 765, groups: simp_plot, 164
##
## Fixed effects:
## Estimate Std. Error t value
## (Intercept) 13.0950 0.3456 37.896
## ap -0.3580 2.4652 -0.145
##
## Correlation of Fixed Effects:
## (Intr)
## ap -0.716Rich_bryo<-lmer(S_plotdate~bryo+(1|simp_plot),vegecomm_bycount)
summary(Rich_bryo)## Linear mixed model fit by REML ['lmerMod']
## Formula: S_plotdate ~ bryo + (1 | simp_plot)
## Data: vegecomm_bycount
##
## REML criterion at convergence: 4333.5
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -2.5734 -0.6798 -0.0879 0.6431 3.5273
##
## Random effects:
## Groups Name Variance Std.Dev.
## simp_plot (Intercept) 5.611 2.369
## Residual 13.679 3.699
## Number of obs: 765, groups: simp_plot, 164
##
## Fixed effects:
## Estimate Std. Error t value
## (Intercept) 13.6941 0.2986 45.862
## bryo -42.5099 12.3301 -3.448
##
## Correlation of Fixed Effects:
## (Intr)
## bryo -0.625Rich_sens_rel<-lmer(S_plotdate~sens_rel+(1|simp_plot),vegecomm_bycount)
summary(Rich_sens_rel)## Linear mixed model fit by REML ['lmerMod']
## Formula: S_plotdate ~ sens_rel + (1 | simp_plot)
## Data: vegecomm_bycount
##
## REML criterion at convergence: 4337.3
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -2.5731 -0.7099 -0.1163 0.6427 3.5445
##
## Random effects:
## Groups Name Variance Std.Dev.
## simp_plot (Intercept) 5.461 2.337
## Residual 13.724 3.705
## Number of obs: 765, groups: simp_plot, 164
##
## Fixed effects:
## Estimate Std. Error t value
## (Intercept) 13.7862 0.3107 44.375
## sens_rel -4.8222 1.3559 -3.556
##
## Correlation of Fixed Effects:
## (Intr)
## sens_rel -0.668Rich_tol_rel<-lmer(S_plotdate~tol_rel+(1|simp_plot),vegecomm_bycount)
summary(Rich_tol_rel)## Linear mixed model fit by REML ['lmerMod']
## Formula: S_plotdate ~ tol_rel + (1 | simp_plot)
## Data: vegecomm_bycount
##
## REML criterion at convergence: 4330.8
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -2.6961 -0.6952 -0.1092 0.6514 3.4089
##
## Random effects:
## Groups Name Variance Std.Dev.
## simp_plot (Intercept) 5.137 2.267
## Residual 13.698 3.701
## Number of obs: 765, groups: simp_plot, 164
##
## Fixed effects:
## Estimate Std. Error t value
## (Intercept) 11.8052 0.3544 33.310
## tol_rel 3.8867 0.8554 4.544
##
## Correlation of Fixed Effects:
## (Intr)
## tol_rel -0.769Rich_fqi<-lmer(S_plotdate~fqi+(1|simp_plot),vegecomm_bycount)
summary(Rich_fqi)## Linear mixed model fit by REML ['lmerMod']
## Formula: S_plotdate ~ fqi + (1 | simp_plot)
## Data: vegecomm_bycount
##
## REML criterion at convergence: 4355.3
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -2.5737 -0.6954 -0.1016 0.6440 3.4662
##
## Random effects:
## Groups Name Variance Std.Dev.
## simp_plot (Intercept) 6.134 2.477
## Residual 13.698 3.701
## Number of obs: 765, groups: simp_plot, 164
##
## Fixed effects:
## Estimate Std. Error t value
## (Intercept) 13.83929 0.70282 19.691
## fqi -0.04354 0.03680 -1.183
##
## Correlation of Fixed Effects:
## (Intr)
## fqi -0.940Rich_inv<-lmer(S_plotdate~inv+(1|simp_plot),vegecomm_bycount)
summary(Rich_inv)## Linear mixed model fit by REML ['lmerMod']
## Formula: S_plotdate ~ inv + (1 | simp_plot)
## Data: vegecomm_bycount
##
## REML criterion at convergence: 4346.2
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -2.5113 -0.6818 -0.0954 0.6528 3.5191
##
## Random effects:
## Groups Name Variance Std.Dev.
## simp_plot (Intercept) 6.115 2.473
## Residual 13.697 3.701
## Number of obs: 765, groups: simp_plot, 164
##
## Fixed effects:
## Estimate Std. Error t value
## (Intercept) 13.0095 0.2425 53.658
## inv 3.7372 2.8028 1.333
##
## Correlation of Fixed Effects:
## (Intr)
## inv -0.148Rich_small_tree<-lmer(S_plotdate~small_tree+(1|simp_plot),vegecomm_bycount)
summary(Rich_small_tree)## Linear mixed model fit by REML ['lmerMod']
## Formula: S_plotdate ~ small_tree + (1 | simp_plot)
## Data: vegecomm_bycount
##
## REML criterion at convergence: 4348.3
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -2.5315 -0.6775 -0.0946 0.6508 3.4892
##
## Random effects:
## Groups Name Variance Std.Dev.
## simp_plot (Intercept) 6.242 2.498
## Residual 13.688 3.700
## Number of obs: 765, groups: simp_plot, 164
##
## Fixed effects:
## Estimate Std. Error t value
## (Intercept) 13.1300 0.4135 31.750
## small_tree -0.4937 2.3432 -0.211
##
## Correlation of Fixed Effects:
## (Intr)
## small_tree -0.812Rich_subcanopy<-lmer(S_plotdate~subcanopy+(1|simp_plot),vegecomm_bycount)
summary(Rich_subcanopy)## Linear mixed model fit by REML ['lmerMod']
## Formula: S_plotdate ~ subcanopy + (1 | simp_plot)
## Data: vegecomm_bycount
##
## REML criterion at convergence: 4348.5
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -2.5203 -0.6845 -0.0962 0.6450 3.4966
##
## Random effects:
## Groups Name Variance Std.Dev.
## simp_plot (Intercept) 6.238 2.498
## Residual 13.689 3.700
## Number of obs: 765, groups: simp_plot, 164
##
## Fixed effects:
## Estimate Std. Error t value
## (Intercept) 13.0349 0.2907 44.85
## subcanopy 0.3244 2.1587 0.15
##
## Correlation of Fixed Effects:
## (Intr)
## subcanopy -0.557Rich_canopy<-lmer(S_plotdate~canopy+(1|simp_plot),vegecomm_bycount)
summary(Rich_canopy)## Linear mixed model fit by REML ['lmerMod']
## Formula: S_plotdate ~ canopy + (1 | simp_plot)
## Data: vegecomm_bycount
##
## REML criterion at convergence: 4340.5
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -2.7009 -0.6938 -0.0975 0.6317 3.6263
##
## Random effects:
## Groups Name Variance Std.Dev.
## simp_plot (Intercept) 5.823 2.413
## Residual 13.686 3.699
## Number of obs: 765, groups: simp_plot, 164
##
## Fixed effects:
## Estimate Std. Error t value
## (Intercept) 14.7055 0.6461 22.762
## canopy -7.9969 2.9116 -2.747
##
## Correlation of Fixed Effects:
## (Intr)
## canopy -0.931Rich_null<-lmer(S_plotdate~1+(1|simp_plot),vegecomm_bycount)
summary(Rich_null)## Linear mixed model fit by REML ['lmerMod']
## Formula: S_plotdate ~ 1 + (1 | simp_plot)
## Data: vegecomm_bycount
##
## REML criterion at convergence: 4351.9
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -2.5202 -0.6845 -0.0981 0.6522 3.5041
##
## Random effects:
## Groups Name Variance Std.Dev.
## simp_plot (Intercept) 6.176 2.485
## Residual 13.691 3.700
## Number of obs: 765, groups: simp_plot, 164
##
## Fixed effects:
## Estimate Std. Error t value
## (Intercept) 13.0584 0.2406 54.28richmods<-list(Rich_ap,Rich_bryo,Rich_canopy,Rich_carex,Rich_cyper,Rich_dicot,Rich_fqi,Rich_hydro,Rich_inv,Rich_sens_rel,Rich_shade,Rich_shrub,Rich_small_tree,Rich_subcanopy,Rich_svp,Rich_tol_rel,Rich_vibifq,Rich_vibifqnotrees,Rich_null)
aictab(richmods)##
## Model selection based on AICc:
##
## K AICc Delta_AICc AICcWt Cum.Wt Res.LL
## Mod16 4 4338.87 0.00 0.76 0.76 -2165.41
## Mod2 4 4341.52 2.65 0.20 0.96 -2166.73
## Mod10 4 4345.35 6.48 0.03 0.99 -2168.65
## Mod3 4 4348.56 9.69 0.01 1.00 -2170.26
## Mod9 4 4354.29 15.43 0.00 1.00 -2173.12
## Mod1 4 4356.31 17.44 0.00 1.00 -2174.13
## Mod13 4 4356.39 17.52 0.00 1.00 -2174.17
## Mod14 4 4356.57 17.70 0.00 1.00 -2174.26
## Mod8 4 4357.41 18.54 0.00 1.00 -2174.68
## Mod19 3 4357.95 19.08 0.00 1.00 -2175.96
## Mod17 4 4358.70 19.83 0.00 1.00 -2175.33
## Mod12 4 4359.36 20.50 0.00 1.00 -2175.66
## Mod6 4 4359.50 20.63 0.00 1.00 -2175.72
## Mod5 4 4361.74 22.87 0.00 1.00 -2176.84
## Mod15 4 4361.86 22.99 0.00 1.00 -2176.90
## Mod4 4 4361.86 22.99 0.00 1.00 -2176.90
## Mod7 4 4363.34 24.47 0.00 1.00 -2177.64
## Mod11 4 4364.24 25.37 0.00 1.00 -2178.09
## Mod18 4 4364.98 26.11 0.00 1.00 -2178.46So richness - at a plot level, and after accounting for plot identity - is determined by how many tolerant species are in the associated plot. And this is much more predictive than any other model.
Here is the plot for this top model:
## `geom_smooth()` using formula 'y ~ x'
We can also verify that our mixed effects model fits better compared to the fixed effects model.
Rich_tol_rel_fixed<-lm(S_plotdate~tol_rel,vegecomm_bycount)
summary(Rich_tol_rel_fixed)##
## Call:
## lm(formula = S_plotdate ~ tol_rel, data = vegecomm_bycount)
##
## Residuals:
## Min 1Q Median 3Q Max
## -10.515 -2.830 -0.563 2.757 19.963
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 11.6653 0.2401 48.585 < 2e-16 ***
## tol_rel 3.7349 0.5971 6.255 6.6e-10 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 4.29 on 763 degrees of freedom
## Multiple R-squared: 0.04878, Adjusted R-squared: 0.04754
## F-statistic: 39.13 on 1 and 763 DF, p-value: 6.604e-10anova(Rich_tol_rel, Rich_tol_rel_fixed) ## refitting model(s) with ML (instead of REML)## Data: vegecomm_bycount
## Models:
## Rich_tol_rel_fixed: S_plotdate ~ tol_rel
## Rich_tol_rel: S_plotdate ~ tol_rel + (1 | simp_plot)
## npar AIC BIC logLik deviance Chisq Df Pr(>Chisq)
## Rich_tol_rel_fixed 3 4403.2 4417.1 -2198.6 4397.2
## Rich_tol_rel 4 4339.2 4357.8 -2165.6 4331.2 66.02 1 4.464e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1This confirms that we chose correctly - the mixed effects model fits better than the fixed effects-only model.
We can also inspect how this observation-level, mixed effect model compares to a more consolidated, plot-level view (in other words - a view that just looks at richness at a plot level - aggregating over observations).
S_tol_rel<-lm(S~tol_rel,vegecomm_richness)
confint(S_tol_rel)## 2.5 % 97.5 %
## (Intercept) 22.889198 25.86972
## tol_rel 1.996856 9.10019confint(Rich_tol_rel)## Computing profile confidence intervals ...## 2.5 % 97.5 %
## .sig01 1.832966 2.688447
## .sigma 3.500278 3.923018
## (Intercept) 11.111432 12.501114
## tol_rel 2.211719 5.563805We can see that, while the covariate estimate for the plot-level model is higher, it is also estimated with much more variance (its confidence intervals are wider). This is not surprising - our functional sample size for this dataset is 164 plots. In contrast, while the observation-level mixed-effects model has a lower covariate estimate, we have tighter confidence intervals around this estimate - likely by virtue of a larger number of samples and a more effectively parameterized model.
reservation-level analyses
birdcomm_reserv<-CMP_bbs_pcap_join_pcaponly %>%
group_by(reservatio,species_4_) %>%
summarise(n=n()) %>%
spread(species_4_,n) %>%
replace(is.na(.),0) ## `summarise()` has grouped output by 'reservatio'. You can override using the `.groups` argument.birdcomm_reserv<-column_to_rownames(birdcomm_reserv,var="reservatio")
birdcomm_reserv<-birdcomm_reserv[,-1]
birdcomm_reserv_bin<-birdcomm_reserv %>% mutate_if(is.numeric, ~1 * (. >= 1))
reservspcount<-specnumber(birdcomm_reserv_bin)
reservspcount## BC BE BK BR BW EC GP HI HU MS NC OE RR SC WA WC
## 48 73 45 70 41 20 24 78 27 63 66 43 80 76 23 48vegecomm_reserv<-CMP_bbs_pcap_join_pcaponly %>%
group_by(reservatio) %>%
summarise(vibifq=mean(vibi_fq),vibifqnotrees=mean(vibi_fq_no_trees),carex=mean(carex_metric_value),cyper=mean(cyperaceae_metric_value),
dicot=mean(dicot_metric_value),shade=mean(shade_metric_value),shrub=mean(shrub_metric_value),
hydro=mean(hydrophyte_metric_value),svp=mean(svp_metric_value),ap=mean(ap_ratio_metric_value),
fqi=mean(fqai_score),bryo=mean(bryophyte_metric_value),hydro_rel=mean(hydrophyte_rel_cov_metric_value),
hydro_rel_notrees=mean(hydrophyte_rel_cov_no_trees),
sens_rel=mean(sensitive_rel_cov_metric_value),
sens_rel_notree=mean(sensitive_rel_cov_no_trees),
tol_rel=mean(tolerant_rel_cov_metric_value),
tol_rel_notree=mean(tolerant_rel_cov_no_trees),
inv=mean(invasive_graminoids_metric_value),
small_tree=mean(small_tree_metric_value),
subcanopy=mean(subcanopy_iv),
canopy=mean(canopy_iv.1),
numsurveys=n_distinct(paste0(simp_plot,date)))
vegecomm_reserv_richness<-cbind(vegecomm_reserv,reservspcount)
cor(vegecomm_reserv_richness[,-1])## vibifq vibifqnotrees carex cyper dicot shade shrub hydro svp ap fqi bryo hydro_rel hydro_rel_notrees sens_rel sens_rel_notree tol_rel tol_rel_notree inv small_tree subcanopy
## vibifq 1.00000000 0.9241914 0.50266385 0.46630888 0.564528005 0.65220135 0.7527538 0.62273545 0.60880098 0.257794968 0.947059034 0.30558082 0.62339581 0.65183940 0.56704683 0.42599572 -0.852896083 -0.8105305 0.32149298 0.47496212 0.531367247
## vibifqnotrees 0.92419139 1.0000000 0.53929085 0.47464154 0.449543250 0.61107030 0.7782038 0.55620231 0.69989878 0.418705980 0.874506913 0.36484997 0.59042482 0.60034564 0.52764476 0.53330180 -0.774869734 -0.8047408 0.28706233 0.40312084 0.498283056
## carex 0.50266385 0.5392908 1.00000000 0.98861968 0.572023490 0.75776669 0.5042649 0.35192082 0.64080449 0.377693834 0.671921171 -0.25412469 0.32776255 0.34946703 0.07325813 -0.08076780 -0.304555329 -0.2892609 0.53268817 0.18937253 -0.281897425
## cyper 0.46630888 0.4746415 0.98861968 1.00000000 0.584049479 0.75372315 0.4711865 0.30703125 0.58053080 0.380730464 0.647650588 -0.27287406 0.31028593 0.32971714 0.06376714 -0.10310014 -0.258492814 -0.2380678 0.50561201 0.17547331 -0.327850311
## dicot 0.56452800 0.4495432 0.57202349 0.58404948 1.000000000 0.94021548 0.4698399 0.54034866 0.48543907 0.008029017 0.759772038 -0.24608019 0.28671560 0.24788748 -0.08333476 -0.19172516 -0.347413659 -0.1282112 0.17164936 -0.03846365 0.155163341
## shade 0.65220135 0.6110703 0.75776669 0.75372315 0.940215480 1.00000000 0.5829774 0.52014187 0.66723418 0.142508957 0.832756899 -0.22197887 0.35664035 0.33764329 -0.01115955 -0.10193730 -0.369664790 -0.2287267 0.35303648 0.06045936 0.084086166
## shrub 0.75275378 0.7782038 0.50426492 0.47118653 0.469839900 0.58297741 1.0000000 0.78513894 0.61806333 0.111257132 0.750745575 -0.13100838 0.83342517 0.83262041 0.23501494 0.14386524 -0.460370277 -0.4995103 0.10638055 0.41413586 0.249525499
## hydro 0.62273545 0.5562023 0.35192082 0.30703125 0.540348664 0.52014187 0.7851389 1.00000000 0.40397499 -0.209004300 0.622184752 -0.17324938 0.62466824 0.59390990 0.20199468 0.01283236 -0.409560618 -0.3488116 0.12178203 0.20959270 0.183422757
## svp 0.60880098 0.6998988 0.64080449 0.58053080 0.485439069 0.66723418 0.6180633 0.40397499 1.00000000 0.124468114 0.688362107 -0.04239540 0.33225582 0.38825400 0.11004921 0.13158407 -0.319130764 -0.3035986 0.33680986 0.06588272 0.216164439
## ap 0.25779497 0.4187060 0.37769383 0.38073046 0.008029017 0.14250896 0.1112571 -0.20900430 0.12446811 1.000000000 0.268493012 0.46101347 0.16767876 0.14567773 0.42055027 0.56317362 -0.425818603 -0.4889413 0.04450250 0.09192957 0.037297463
## fqi 0.94705903 0.8745069 0.67192117 0.64765059 0.759772038 0.83275690 0.7507456 0.62218475 0.68836211 0.268493012 1.000000000 0.12111170 0.57120898 0.58946433 0.39402289 0.24845247 -0.761357994 -0.6636032 0.30819955 0.33480585 0.391252116
## bryo 0.30558082 0.3648500 -0.25412469 -0.27287406 -0.246080192 -0.22197887 -0.1310084 -0.17324938 -0.04239540 0.461013466 0.121111704 1.00000000 -0.11350270 -0.10247675 0.76072190 0.92221351 -0.529298731 -0.5976709 -0.02147223 0.02772700 0.472984957
## hydro_rel 0.62339581 0.5904248 0.32776255 0.31028593 0.286715602 0.35664035 0.8334252 0.62466824 0.33225582 0.167678756 0.571208983 -0.11350270 1.00000000 0.97600979 0.34652164 0.09194741 -0.511414589 -0.5937322 -0.10008676 0.61633226 0.193693179
## hydro_rel_notrees 0.65183940 0.6003456 0.34946703 0.32971714 0.247887485 0.33764329 0.8326204 0.59390990 0.38825400 0.145677729 0.589464327 -0.10247675 0.97600979 1.00000000 0.39250298 0.10495976 -0.534687160 -0.6144103 -0.01962444 0.58390153 0.178383496
## sens_rel 0.56704683 0.5276448 0.07325813 0.06376714 -0.083334765 -0.01115955 0.2350149 0.20199468 0.11004921 0.420550273 0.394022888 0.76072190 0.34652164 0.39250298 1.00000000 0.82206470 -0.684201135 -0.7837929 0.06059223 0.21381404 0.222296394
## sens_rel_notree 0.42599572 0.5333018 -0.08076780 -0.10310014 -0.191725160 -0.10193730 0.1438652 0.01283236 0.13158407 0.563173621 0.248452466 0.92221351 0.09194741 0.10495976 0.82206470 1.00000000 -0.518781157 -0.6364926 0.01042068 0.05071249 0.383357352
## tol_rel -0.85289608 -0.7748697 -0.30455533 -0.25849281 -0.347413659 -0.36966479 -0.4603703 -0.40956062 -0.31913076 -0.425818603 -0.761357994 -0.52929873 -0.51141459 -0.53468716 -0.68420113 -0.51878116 1.000000000 0.9366517 -0.18307775 -0.45979912 -0.580636269
## tol_rel_notree -0.81053054 -0.8047408 -0.28926090 -0.23806784 -0.128211218 -0.22872670 -0.4995103 -0.34881162 -0.30359859 -0.488941346 -0.663603226 -0.59767089 -0.59373217 -0.61441028 -0.78379286 -0.63649263 0.936651732 1.0000000 -0.16840035 -0.58104234 -0.502322014
## inv 0.32149298 0.2870623 0.53268817 0.50561201 0.171649356 0.35303648 0.1063806 0.12178203 0.33680986 0.044502499 0.308199552 -0.02147223 -0.10008676 -0.01962444 0.06059223 0.01042068 -0.183077752 -0.1684003 1.00000000 0.23345191 -0.168146221
## small_tree 0.47496212 0.4031208 0.18937253 0.17547331 -0.038463654 0.06045936 0.4141359 0.20959270 0.06588272 0.091929570 0.334805846 0.02772700 0.61633226 0.58390153 0.21381404 0.05071249 -0.459799124 -0.5810423 0.23345191 1.00000000 0.306801643
## subcanopy 0.53136725 0.4982831 -0.28189742 -0.32785031 0.155163341 0.08408617 0.2495255 0.18342276 0.21616444 0.037297463 0.391252116 0.47298496 0.19369318 0.17838350 0.22229639 0.38335735 -0.580636269 -0.5023220 -0.16814622 0.30680164 1.000000000
## canopy 0.02686002 -0.1477253 -0.16296479 -0.12875828 0.004174685 -0.07217575 0.1982734 0.12798031 0.01062822 -0.321878438 -0.009884374 -0.27391776 0.52889987 0.54700481 0.01841165 -0.24518828 0.006449976 -0.0197229 -0.24910337 0.35704003 -0.006191979
## numsurveys 0.67749638 0.7157262 0.43821784 0.38388893 0.511769451 0.64351866 0.7158910 0.51650992 0.66735153 0.213173431 0.671919649 -0.06541184 0.60459521 0.63984087 0.15097169 0.10308858 -0.480145398 -0.4446043 0.30324470 0.21309096 0.354918723
## reservspcount 0.68233393 0.7884269 0.50420389 0.45376839 0.464086302 0.60821512 0.7523661 0.54918263 0.62534630 0.449063177 0.673416476 0.08143411 0.66313940 0.63780487 0.24936734 0.27304031 -0.548039099 -0.5601424 0.23875891 0.30577390 0.344432018
## canopy numsurveys reservspcount
## vibifq 0.026860021 0.67749638 0.68233393
## vibifqnotrees -0.147725299 0.71572622 0.78842694
## carex -0.162964792 0.43821784 0.50420389
## cyper -0.128758278 0.38388893 0.45376839
## dicot 0.004174685 0.51176945 0.46408630
## shade -0.072175750 0.64351866 0.60821512
## shrub 0.198273381 0.71589097 0.75236610
## hydro 0.127980314 0.51650992 0.54918263
## svp 0.010628216 0.66735153 0.62534630
## ap -0.321878438 0.21317343 0.44906318
## fqi -0.009884374 0.67191965 0.67341648
## bryo -0.273917763 -0.06541184 0.08143411
## hydro_rel 0.528899874 0.60459521 0.66313940
## hydro_rel_notrees 0.547004812 0.63984087 0.63780487
## sens_rel 0.018411651 0.15097169 0.24936734
## sens_rel_notree -0.245188281 0.10308858 0.27304031
## tol_rel 0.006449976 -0.48014540 -0.54803910
## tol_rel_notree -0.019722898 -0.44460427 -0.56014242
## inv -0.249103369 0.30324470 0.23875891
## small_tree 0.357040034 0.21309096 0.30577390
## subcanopy -0.006191979 0.35491872 0.34443202
## canopy 1.000000000 0.07617249 -0.01191566
## numsurveys 0.076172489 1.00000000 0.90376963
## reservspcount -0.011915663 0.90376963 1.00000000ggplot(aes(x=vibifq,y=reservspcount),data=vegecomm_reserv_richness)+geom_point()+geom_smooth(method="lm")## `geom_smooth()` using formula 'y ~ x'
reservation_spaccum<-specaccum(birdcomm_reserv_bin,w = vegecomm_reserv$numsurveys, method="random")
plot(reservation_spaccum)
This doesn’t quite get at what we want - this displays species accumulation across reservations after accounting for survey effort per reservation. Let’s try exploring reservation-by-reservation rarefaction with the previous observation-level dataset.
birds_bycount## # A tibble: 765 x 104
## ACFL ALFL AMCR AMGO AMKE AMRE AMRO AMWO BADO BAEA BANS BAOR BARS BCCH BEKI BGGN BHCO BHVI BLJA BOBO BRCR BRTH BTNW BWHA BWWA CANG CARW CEDW CERW CHSP CHSW CLSW COGR COHA COYE DEJU DOWO EABL EAKI EAME EAPH EATO EAWP EUST FISP GBHE GCFL GHOW GRCA
## <int> <int> <int> <int> <int> <int> <int> <int> <int> <int> <int> <int> <int> <int> <int> <int> <int> <int> <int> <int> <int> <int> <int> <int> <int> <int> <int> <int> <int> <int> <int> <int> <int> <int> <int> <int> <int> <int> <int> <int> <int> <int> <int> <int> <int> <int> <int> <int> <int>
## 1 0 0 1 1 0 1 0 0 0 0 0 0 0 0 0 1 0 0 1 0 0 0 0 0 1 0 0 1 0 0 0 0 0 0 0 0 1 0 0 0 0 1 1 0 0 0 0 0 0
## 2 0 0 0 1 0 0 0 0 0 0 0 0 0 0 0 0 1 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 1 0 0 0 0 1 1 0 0 0 0 0 0
## 3 1 0 0 0 0 1 0 0 0 0 0 0 0 1 0 0 1 0 0 0 0 0 0 0 0 0 0 0 0 0 1 0 0 0 0 0 0 0 0 0 0 1 0 0 0 0 0 0 0
## 4 0 0 0 1 0 0 1 0 0 0 0 1 0 1 0 0 1 0 1 0 0 0 0 0 0 0 0 1 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0
## 5 1 0 0 1 0 0 1 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 1 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0
## 6 1 0 0 1 0 0 0 0 0 0 0 0 0 0 0 0 1 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 1 0 0
## 7 1 0 0 0 0 0 1 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 1 0 0 0 0 0 0
## 8 1 0 0 0 0 0 0 0 0 0 0 1 0 1 0 1 0 0 1 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 1 0 0 0 0 0 1 0 0 0 0 0 0
## 9 1 0 0 1 0 0 1 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 1 0 0 0 0 0 0 0 0 0 1 0 0 0 0 0 0
## 10 1 0 1 0 0 0 1 0 0 0 0 0 0 1 1 0 0 0 1 0 0 0 0 0 0 0 0 0 0 0 1 0 0 0 0 0 1 0 0 0 0 0 1 0 0 0 1 0 0
## # ... with 755 more rows, and 55 more variables: GRHE <int>, HAWO <int>, HERG <int>, HETH <int>, HOFI <int>, HOSP <int>, HOWA <int>, HOWR <int>, INBU <int>, KILL <int>, LOWA <int>, MALL <int>, MODO <int>, NOCA <int>, NOFL <int>, NOPA <int>, NRWS <int>, OROR <int>, OVEN <int>, PIWO <int>,
## # PUFI <int>, PUMA <int>, RBGR <int>, RBGU <int>, RBNU <int>, RBWO <int>, REVI <int>, RHWO <int>, ROPI <int>, RSHA <int>, RTHA <int>, RTHU <int>, RWBL <int>, SAVS <int>, SCTA <int>, SOSP <int>, SWSP <int>, TEWA <int>, TRES <int>, TUTI <int>, TUVU <int>, VEER <int>, VIRA <int>, WAVI <int>,
## # WBNU <int>, WEVI <int>, WIFL <int>, WITU <int>, WIWR <int>, WODU <int>, WOTH <int>, YBCU <int>, YEWA <int>, YTVI <int>, YTWA <int>vegecomm_bycount## # A tibble: 765 x 28
## # Groups: simp_plot, reservatio, date [765]
## simp_plot reservatio date longitude latitude vibifq vibifqnotrees carex cyper dicot shade shrub hydro svp ap fqi bryo hydro_rel hydro_rel_notrees sens_rel sens_rel_notree tol_rel tol_rel_notree inv small_tree subcanopy canopy V1
## <int> <chr> <chr> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <int>
## 1 1002 HI 6/18/2019 2181968. 561931. 47.3 24.4 9 9 30 25 1 9 1 0.0588 19.5 0.00527 0.0158 0.0317 0.0396 0.000314 0.419 0.766 0 0.139 0.0342 0.192 0
## 2 1002 HI 6/25/2019 2181968. 561931. 47.3 24.4 9 9 30 25 1 9 1 0.0588 19.5 0.00527 0.0158 0.0317 0.0396 0.000314 0.419 0.766 0 0.139 0.0342 0.192 0
## 3 1002 HI 6/5/2019 2181968. 561931. 47.3 24.4 9 9 30 25 1 9 1 0.0588 19.5 0.00527 0.0158 0.0317 0.0396 0.000314 0.419 0.766 0 0.139 0.0342 0.192 0
## 4 1003 MS 6/13/2017 2150599. 613095. 76.2 52.7 4 4 22 22 1 2 1 0.133 23.0 0.00903 0.0181 0.0668 0.488 0.0517 0.111 0.0428 0 0.0485 0.0350 0.205 0
## 5 1003 MS 6/24/2018 2150599. 613095. 76.2 52.7 4 4 22 22 1 2 1 0.133 23.0 0.00903 0.0181 0.0668 0.488 0.0517 0.111 0.0428 0 0.0485 0.0350 0.205 0
## 6 1003 MS 6/25/2017 2150599. 613095. 76.2 52.7 4 4 22 22 1 2 1 0.133 23.0 0.00903 0.0181 0.0668 0.488 0.0517 0.111 0.0428 0 0.0485 0.0350 0.205 0
## 7 1003 MS 6/4/2017 2150599. 613095. 76.2 52.7 4 4 22 22 1 2 1 0.133 23.0 0.00903 0.0181 0.0668 0.488 0.0517 0.111 0.0428 0 0.0485 0.0350 0.205 1
## 8 1003 MS 6/7/2018 2150599. 613095. 76.2 52.7 4 4 22 22 1 2 1 0.133 23.0 0.00903 0.0181 0.0668 0.488 0.0517 0.111 0.0428 0 0.0485 0.0350 0.205 0
## 9 1003 MS 7/1/2018 2150599. 613095. 76.2 52.7 4 4 22 22 1 2 1 0.133 23.0 0.00903 0.0181 0.0668 0.488 0.0517 0.111 0.0428 0 0.0485 0.0350 0.205 0
## 10 1004 BE 6/12/2018 2227458. 626184. 37.9 35.1 2 2 21 23 1 9 4 0.190 15.1 0.0396 0.426 0.559 0.0140 0.00916 0.199 0.260 0 0.109 0.464 0.217 0
## # ... with 755 more rowsindex<-as.numeric(factor(CMP_bbs_pcap_join_pcaponly_pcapsum_date$reservatio))
index_unique<-unique(index)
rarefylist_reservation<-list()
for(i in index_unique){
birds_temp<-birds_bycount[index %in% i,]
rarefylist_reservation[[i]]<- specaccum(birds_temp)
}## Warning in cor(x > 0): the standard deviation is zero
## Warning in cor(x > 0): the standard deviation is zero
## Warning in cor(x > 0): the standard deviation is zero
## Warning in cor(x > 0): the standard deviation is zero
## Warning in cor(x > 0): the standard deviation is zeroplot(totalspecaccum)
plot(rarefylist_reservation[[4]])
interpolatedrichness_n3<-data.frame(levels(factor(CMP_bbs_pcap_join_pcaponly_pcapsum_date$reservatio)),rep(1,16))
colnames(interpolatedrichness_n3)<-c("reservation","richness")
for(i in index_unique){
birds_temp<-birds_bycount[index %in% i,]
interpolatedrichness_n3$richness[i]<- predict(specaccum(birds_temp),newdata=3)
}## Warning in cor(x > 0): the standard deviation is zero
## Warning in cor(x > 0): the standard deviation is zero
## Warning in cor(x > 0): the standard deviation is zero
## Warning in cor(x > 0): the standard deviation is zero
## Warning in cor(x > 0): the standard deviation is zeroWith the plot-level data, we could easily account for variation in survey effort across plots with a mixed model design, but this case is a bit more complicated. We are interested in reservation-level species richness, but know that reservations were visited different numbers of times - based on the number of plots and the number of surveys/plot. Thus, the easiest thing to start with is simply using species richness estimates interpolated to the least number of surveys per reservation (3, in our case).
Now we can proceed with exploring richness patterns across reservations.
vegecomm_reserv_richness$interrich_n3<-interpolatedrichness_n3$richness
Rich_reserv_vibifq<-lm(interrich_n3~vibifq,vegecomm_reserv_richness)
summary(Rich_reserv_vibifq)##
## Call:
## lm(formula = interrich_n3 ~ vibifq, data = vegecomm_reserv_richness)
##
## Residuals:
## Min 1Q Median 3Q Max
## -5.0917 -1.9946 0.4847 1.9275 4.8225
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 19.24757 1.87384 10.272 6.69e-08 ***
## vibifq 0.12874 0.04718 2.729 0.0163 *
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 2.731 on 14 degrees of freedom
## Multiple R-squared: 0.3472, Adjusted R-squared: 0.3005
## F-statistic: 7.445 on 1 and 14 DF, p-value: 0.01632Rich_reserv_vibifqnotrees<-lm(interrich_n3~vibifqnotrees,vegecomm_reserv_richness)
summary(Rich_reserv_vibifqnotrees)##
## Call:
## lm(formula = interrich_n3 ~ vibifqnotrees, data = vegecomm_reserv_richness)
##
## Residuals:
## Min 1Q Median 3Q Max
## -5.3776 -1.4771 0.0605 1.4821 4.7685
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 19.81963 1.40534 14.103 1.15e-09 ***
## vibifqnotrees 0.16551 0.04961 3.336 0.0049 **
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 2.523 on 14 degrees of freedom
## Multiple R-squared: 0.4429, Adjusted R-squared: 0.4031
## F-statistic: 11.13 on 1 and 14 DF, p-value: 0.004899Rich_reserv_carex<-lm(interrich_n3~carex,vegecomm_reserv_richness)
summary(Rich_reserv_carex)##
## Call:
## lm(formula = interrich_n3 ~ carex, data = vegecomm_reserv_richness)
##
## Residuals:
## Min 1Q Median 3Q Max
## -4.7041 -2.3742 0.1555 1.4905 5.9512
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 22.2721 1.6801 13.257 2.58e-09 ***
## carex 0.6488 0.5506 1.178 0.258
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 3.224 on 14 degrees of freedom
## Multiple R-squared: 0.09023, Adjusted R-squared: 0.02525
## F-statistic: 1.388 on 1 and 14 DF, p-value: 0.2583Rich_reserv_cyper<-lm(interrich_n3~cyper,vegecomm_reserv_richness)
summary(Rich_reserv_cyper)##
## Call:
## lm(formula = interrich_n3 ~ cyper, data = vegecomm_reserv_richness)
##
## Residuals:
## Min 1Q Median 3Q Max
## -5.091 -2.482 0.218 1.558 5.838
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 22.5436 1.7022 13.243 2.61e-09 ***
## cyper 0.5155 0.5252 0.982 0.343
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 3.269 on 14 degrees of freedom
## Multiple R-squared: 0.06438, Adjusted R-squared: -0.00245
## F-statistic: 0.9633 on 1 and 14 DF, p-value: 0.343Rich_reserv_dicot<-lm(interrich_n3~dicot,vegecomm_reserv_richness)
summary(Rich_reserv_dicot)##
## Call:
## lm(formula = interrich_n3 ~ dicot, data = vegecomm_reserv_richness)
##
## Residuals:
## Min 1Q Median 3Q Max
## -5.896 -1.993 0.058 1.958 4.709
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 20.9424 2.2412 9.344 2.15e-07 ***
## dicot 0.1489 0.1019 1.461 0.166
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 3.148 on 14 degrees of freedom
## Multiple R-squared: 0.1324, Adjusted R-squared: 0.07039
## F-statistic: 2.136 on 1 and 14 DF, p-value: 0.166Rich_reserv_shade<-lm(interrich_n3~shade,vegecomm_reserv_richness)
summary(Rich_reserv_shade)##
## Call:
## lm(formula = interrich_n3 ~ shade, data = vegecomm_reserv_richness)
##
## Residuals:
## Min 1Q Median 3Q Max
## -5.532 -1.853 -0.095 1.742 4.954
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 21.09296 1.91432 11.018 2.77e-08 ***
## shade 0.16474 0.09896 1.665 0.118
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 3.088 on 14 degrees of freedom
## Multiple R-squared: 0.1652, Adjusted R-squared: 0.1056
## F-statistic: 2.771 on 1 and 14 DF, p-value: 0.1182Rich_reserv_hydro<-lm(interrich_n3~hydro,vegecomm_reserv_richness)
summary(Rich_reserv_hydro)##
## Call:
## lm(formula = interrich_n3 ~ hydro, data = vegecomm_reserv_richness)
##
## Residuals:
## Min 1Q Median 3Q Max
## -3.3095 -2.0680 0.4451 1.0590 5.1184
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 21.3598 1.1362 18.800 2.48e-11 ***
## hydro 0.6499 0.2252 2.885 0.012 *
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 2.676 on 14 degrees of freedom
## Multiple R-squared: 0.3729, Adjusted R-squared: 0.3281
## F-statistic: 8.324 on 1 and 14 DF, p-value: 0.01199Rich_reserv_shrub<-lm(interrich_n3~shrub,vegecomm_reserv_richness)
summary(Rich_reserv_shrub)##
## Call:
## lm(formula = interrich_n3 ~ shrub, data = vegecomm_reserv_richness)
##
## Residuals:
## Min 1Q Median 3Q Max
## -4.5855 -2.2688 0.3838 1.5842 4.6146
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 22.147 1.055 20.996 5.56e-12 ***
## shrub 6.154 2.573 2.392 0.0314 *
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 2.848 on 14 degrees of freedom
## Multiple R-squared: 0.2901, Adjusted R-squared: 0.2394
## F-statistic: 5.721 on 1 and 14 DF, p-value: 0.03135Rich_reserv_svp<-lm(interrich_n3~svp,vegecomm_reserv_richness)
summary(Rich_reserv_svp)##
## Call:
## lm(formula = interrich_n3 ~ svp, data = vegecomm_reserv_richness)
##
## Residuals:
## Min 1Q Median 3Q Max
## -4.8940 -2.3474 0.0939 2.4361 5.6107
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 23.2440 1.1032 21.069 5.3e-12 ***
## svp 0.9503 0.9240 1.028 0.321
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 3.259 on 14 degrees of freedom
## Multiple R-squared: 0.07025, Adjusted R-squared: 0.003838
## F-statistic: 1.058 on 1 and 14 DF, p-value: 0.3212Rich_reserv_ap<-lm(interrich_n3~ap,vegecomm_reserv_richness)
summary(Rich_reserv_ap)##
## Call:
## lm(formula = interrich_n3 ~ ap, data = vegecomm_reserv_richness)
##
## Residuals:
## Min 1Q Median 3Q Max
## -4.8181 -2.0705 -0.2523 1.5602 4.7442
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 21.675 1.532 14.152 1.1e-09 ***
## ap 25.248 14.352 1.759 0.1
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 3.058 on 14 degrees of freedom
## Multiple R-squared: 0.181, Adjusted R-squared: 0.1225
## F-statistic: 3.095 on 1 and 14 DF, p-value: 0.1004Rich_reserv_bryo<-lm(interrich_n3~bryo,vegecomm_reserv_richness)
summary(Rich_reserv_bryo)##
## Call:
## lm(formula = interrich_n3 ~ bryo, data = vegecomm_reserv_richness)
##
## Residuals:
## Min 1Q Median 3Q Max
## -4.7465 -1.9194 0.1346 1.4286 5.1210
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 23.096 1.047 22.065 2.83e-12 ***
## bryo 59.466 44.336 1.341 0.201
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 3.181 on 14 degrees of freedom
## Multiple R-squared: 0.1139, Adjusted R-squared: 0.05057
## F-statistic: 1.799 on 1 and 14 DF, p-value: 0.2012Rich_reserv_sens_rel<-lm(interrich_n3~sens_rel,vegecomm_reserv_richness)
summary(Rich_reserv_sens_rel)##
## Call:
## lm(formula = interrich_n3 ~ sens_rel, data = vegecomm_reserv_richness)
##
## Residuals:
## Min 1Q Median 3Q Max
## -4.0560 -1.6362 -0.5409 1.4794 5.9762
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 22.042 1.382 15.953 2.25e-10 ***
## sens_rel 14.253 8.320 1.713 0.109
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 3.073 on 14 degrees of freedom
## Multiple R-squared: 0.1733, Adjusted R-squared: 0.1143
## F-statistic: 2.935 on 1 and 14 DF, p-value: 0.1087Rich_reserv_tol_rel<-lm(interrich_n3~tol_rel,vegecomm_reserv_richness)
summary(Rich_reserv_tol_rel)##
## Call:
## lm(formula = interrich_n3 ~ tol_rel, data = vegecomm_reserv_richness)
##
## Residuals:
## Min 1Q Median 3Q Max
## -4.981 -1.519 0.258 1.307 3.891
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 28.411 1.599 17.765 5.33e-11 ***
## tol_rel -10.939 3.623 -3.019 0.00919 **
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 2.63 on 14 degrees of freedom
## Multiple R-squared: 0.3943, Adjusted R-squared: 0.3511
## F-statistic: 9.116 on 1 and 14 DF, p-value: 0.009195Rich_reserv_fqi<-lm(interrich_n3~fqi,vegecomm_reserv_richness)
summary(Rich_reserv_fqi)##
## Call:
## lm(formula = interrich_n3 ~ fqi, data = vegecomm_reserv_richness)
##
## Residuals:
## Min 1Q Median 3Q Max
## -5.1569 -1.5216 0.2207 1.6862 5.0045
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 18.5248 2.4319 7.617 2.41e-06 ***
## fqi 0.3606 0.1528 2.359 0.0334 *
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 2.859 on 14 degrees of freedom
## Multiple R-squared: 0.2845, Adjusted R-squared: 0.2334
## F-statistic: 5.566 on 1 and 14 DF, p-value: 0.03336Rich_reserv_inv<-lm(interrich_n3~inv,vegecomm_reserv_richness)
summary(Rich_reserv_inv)##
## Call:
## lm(formula = interrich_n3 ~ inv, data = vegecomm_reserv_richness)
##
## Residuals:
## Min 1Q Median 3Q Max
## -5.4079 -2.2428 0.2233 1.6034 5.3265
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 23.758 0.951 24.983 5.17e-13 ***
## inv 27.430 49.518 0.554 0.588
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 3.343 on 14 degrees of freedom
## Multiple R-squared: 0.02145, Adjusted R-squared: -0.04845
## F-statistic: 0.3069 on 1 and 14 DF, p-value: 0.5884Rich_reserv_small_tree<-lm(interrich_n3~small_tree,vegecomm_reserv_richness)
summary(Rich_reserv_small_tree)##
## Call:
## lm(formula = interrich_n3 ~ small_tree, data = vegecomm_reserv_richness)
##
## Residuals:
## Min 1Q Median 3Q Max
## -4.9406 -2.1955 0.3979 1.8062 4.2438
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 21.90 1.93 11.346 1.91e-08 ***
## small_tree 15.98 13.30 1.201 0.25
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 3.218 on 14 degrees of freedom
## Multiple R-squared: 0.09346, Adjusted R-squared: 0.02871
## F-statistic: 1.443 on 1 and 14 DF, p-value: 0.2495Rich_reserv_subcanopy<-lm(interrich_n3~subcanopy,vegecomm_reserv_richness)
summary(Rich_reserv_subcanopy)##
## Call:
## lm(formula = interrich_n3 ~ subcanopy, data = vegecomm_reserv_richness)
##
## Residuals:
## Min 1Q Median 3Q Max
## -6.1575 -1.2786 0.9301 1.9275 3.7641
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 22.639 1.269 17.838 5.04e-11 ***
## subcanopy 23.877 17.271 1.383 0.188
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 3.17 on 14 degrees of freedom
## Multiple R-squared: 0.1201, Adjusted R-squared: 0.05727
## F-statistic: 1.911 on 1 and 14 DF, p-value: 0.1885Rich_reserv_canopy<-lm(interrich_n3~canopy,vegecomm_reserv_richness)
summary(Rich_reserv_canopy)##
## Call:
## lm(formula = interrich_n3 ~ canopy, data = vegecomm_reserv_richness)
##
## Residuals:
## Min 1Q Median 3Q Max
## -5.4648 -2.6753 0.5155 1.7726 5.3828
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 26.576 4.432 5.997 3.27e-05 ***
## canopy -12.717 21.561 -0.590 0.565
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 3.338 on 14 degrees of freedom
## Multiple R-squared: 0.02424, Adjusted R-squared: -0.04545
## F-statistic: 0.3479 on 1 and 14 DF, p-value: 0.5647Rich_reserv_null<-lm(interrich_n3~1,vegecomm_reserv_richness)
summary(Rich_reserv_null)##
## Call:
## lm(formula = interrich_n3 ~ 1, data = vegecomm_reserv_richness)
##
## Residuals:
## Min 1Q Median 3Q Max
## -5.6592 -2.4610 0.2916 1.9212 5.1310
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 24.0092 0.8163 29.41 1.11e-14 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 3.265 on 15 degrees of freedomrichmods_reserv<-list(Rich_reserv_ap,Rich_reserv_bryo,Rich_reserv_canopy,Rich_reserv_carex,Rich_reserv_cyper,Rich_reserv_dicot,Rich_reserv_fqi,Rich_reserv_hydro,Rich_reserv_inv,Rich_reserv_sens_rel,Rich_reserv_shade,Rich_reserv_shrub,Rich_reserv_small_tree,Rich_reserv_subcanopy,Rich_reserv_svp,Rich_reserv_tol_rel,Rich_reserv_vibifq,Rich_reserv_vibifqnotrees,Rich_reserv_null)
aictab(richmods_reserv)## Warning in aictab.AIClm(richmods_reserv):
## Model names have been supplied automatically in the table##
## Model selection based on AICc:
##
## K AICc Delta_AICc AICcWt Cum.Wt LL
## Mod18 3 80.88 0.00 0.36 0.36 -36.44
## Mod16 3 82.21 1.34 0.18 0.54 -37.11
## Mod8 3 82.77 1.89 0.14 0.68 -37.39
## Mod17 3 83.41 2.54 0.10 0.78 -37.71
## Mod12 3 84.76 3.88 0.05 0.83 -38.38
## Mod7 3 84.88 4.00 0.05 0.88 -38.44
## Mod1 3 87.04 6.16 0.02 0.89 -39.52
## Mod19 2 87.16 6.28 0.02 0.91 -41.12
## Mod10 3 87.19 6.31 0.02 0.92 -39.60
## Mod11 3 87.35 6.47 0.01 0.94 -39.67
## Mod6 3 87.97 7.09 0.01 0.95 -39.98
## Mod14 3 88.19 7.31 0.01 0.96 -40.10
## Mod2 3 88.30 7.42 0.01 0.97 -40.15
## Mod13 3 88.67 7.79 0.01 0.97 -40.33
## Mod4 3 88.73 7.85 0.01 0.98 -40.36
## Mod15 3 89.07 8.19 0.01 0.99 -40.54
## Mod5 3 89.17 8.29 0.01 0.99 -40.59
## Mod3 3 89.85 8.97 0.00 1.00 -40.92
## Mod9 3 89.89 9.01 0.00 1.00 -40.95although our power is a bit limited above, it does appear that the Vegetative Index of Biotic Integrity - at a reservation-level - does appear to predict bird species richness (as long as we do not include woody plants in that list!)
we can also plot out richness for each site:
ggplot(aes(x=reservatio,y=reservspcount),data=vegecomm_reserv_richness)+geom_col()
ggplot(aes(x=reservatio,y=interrich_n3),data=vegecomm_reserv_richness)+geom_col()
In confirmation of our prior ordination plots, it now appears that Rocky River is one of the most species rich areas - followed by Hinckley, South Chagrin, Bedford, and Brecksville.
reservation shape features and species richness
While these vegetative associations are certainly important, it is often helpful to also think about how the shape of our reservations influences the species found within them. There is a rich body of literature on the effects of landscape context and configuration on
vegecomm_reserv_richness## reservatio vibifq vibifqnotrees carex cyper dicot shade shrub hydro svp ap fqi bryo hydro_rel hydro_rel_notrees sens_rel sens_rel_notree tol_rel tol_rel_notree inv small_tree subcanopy canopy numsurveys reservspcount
## BC BC 34.601035 25.313664 3.8158996 4.1087866 16.669456 15.774059 0.3807531 5.6276151 0.0000000 0.11874312 14.284368 0.004958120 0.023797808 0.028223440 0.089812897 0.045449261 0.4891651 0.5660586 0.0295276552 0.17029270 0.01723615 0.1166671 18 48
## BE BE 39.504613 29.726028 3.0626566 3.1102757 21.440267 20.778613 0.2063492 2.8671679 1.4152047 0.14248544 15.672360 0.014829724 0.047829388 0.078499350 0.164717698 0.070814554 0.3887067 0.4714974 0.0051759616 0.11508660 0.07952218 0.1801468 100 73
## BK BK 40.504890 33.317087 0.6134454 0.6134454 9.773109 6.680672 0.0000000 0.7731092 0.0000000 0.18309632 13.396361 0.081362435 0.000000000 0.000000000 0.366516404 0.310637023 0.1836977 0.1760264 0.0000000000 0.09542484 0.11412427 0.1443458 9 45
## BR BR 55.285843 41.682545 5.2304060 5.4381492 28.100094 30.229462 0.4277620 4.4711992 2.7648725 0.10423918 21.496159 0.022946191 0.035888066 0.070085339 0.166597534 0.103523925 0.3088704 0.3869337 0.0564445807 0.18449116 0.05321523 0.2079372 82 70
## BW BW 47.464208 35.190650 3.4246285 3.5095541 18.643312 18.087049 0.4394904 2.7940552 1.9341826 0.07325213 19.131091 0.016114563 0.049254626 0.106742227 0.160002194 0.044342921 0.2891663 0.3490757 0.0000000000 0.15571638 0.08904453 0.1886476 41 41
## EC EC 31.350954 11.944005 0.0000000 0.0000000 19.000000 10.000000 0.0000000 3.0000000 0.0000000 0.00000000 12.247449 0.012371644 0.000000000 0.000000000 0.024743288 0.000000000 0.4353994 0.6448408 0.0000000000 0.16216216 0.14737364 0.2303164 3 20
## GP GP 27.946924 10.651350 2.0000000 2.0000000 23.000000 16.000000 0.0000000 7.0000000 0.0000000 0.00000000 12.809880 0.009602612 0.016260423 0.025736505 0.177936399 0.000000000 0.4184018 0.6011348 0.0096026120 0.06153846 0.00000000 0.2061372 3 24
## HI HI 56.777233 46.984497 4.2806748 4.3550613 38.047546 34.980061 0.7730061 9.4685583 1.6541411 0.08149579 23.192002 0.009410177 0.088963035 0.127416913 0.107577765 0.055443524 0.2894462 0.3644430 0.0003946167 0.14618568 0.09429889 0.1791128 92 78
## HU HU 24.435647 6.752597 3.0000000 4.0000000 30.000000 24.000000 0.0000000 0.0000000 0.0000000 0.10526316 14.293825 0.004343231 0.000000000 0.000000000 0.078178159 0.000000000 0.6169125 0.8298138 0.0000000000 0.04727273 0.00000000 0.2077735 7 27
## MS MS 46.355466 29.687439 3.5869738 3.7172359 20.918983 19.845909 0.3494837 4.3288324 0.7648928 0.08450301 17.252813 0.011838657 0.048039227 0.123442813 0.165213079 0.032381386 0.2526545 0.3385182 0.0414410260 0.14354510 0.05094113 0.1931117 104 63
## NC NC 53.965315 38.560833 2.3760684 2.4658120 23.248718 20.439316 0.8153846 7.8316239 1.3606838 0.10368316 19.929270 0.014654907 0.112169719 0.209235627 0.227754177 0.127802382 0.2922618 0.3813379 0.0002100554 0.12328842 0.06749561 0.2470918 96 66
## OE OE 13.813003 14.088671 3.6745283 3.6745283 17.641509 14.811321 0.0000000 0.1839623 0.8207547 0.16659520 10.636789 0.002949977 0.001007902 0.001672384 0.005627719 0.006567877 0.5388112 0.6840616 0.0016111062 0.05005438 0.00000000 0.1768963 18 43
## RR RR 39.382635 27.503422 1.4131868 1.4717949 23.423443 18.774359 0.3897436 4.0959707 0.3003663 0.11351880 15.560389 0.015515483 0.084112185 0.123459302 0.078730158 0.025448674 0.2574988 0.3879014 0.0020324549 0.19507519 0.11464307 0.2216719 93 80
## SC SC 39.740571 28.655854 3.1097209 3.3666987 22.470645 18.662175 0.6496631 8.8960539 1.8671800 0.07591598 16.413909 0.010576548 0.072576642 0.128944299 0.147712040 0.073316453 0.4357619 0.5385634 0.0001414818 0.08700562 0.06758911 0.2247417 76 76
## WA WA 1.137611 1.026171 0.0000000 0.0000000 4.000000 3.000000 0.0000000 0.0000000 0.0000000 0.00000000 2.474874 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.9516582 0.9604965 0.0000000000 0.08695652 0.00000000 0.2171573 6 23
## WC WC 39.518844 23.929770 3.2511628 3.6604651 13.069767 11.162791 0.4093023 3.8883721 0.0000000 0.12656059 14.560841 0.014097267 0.126317347 0.209907312 0.247360772 0.065370632 0.2894291 0.2196096 0.0000000000 0.28662562 0.02284111 0.2881646 17 48
## interrich_n3
## BC 25.65319
## BE 25.03640
## BK 26.76190
## BR 25.61177
## BW 20.26632
## EC 20.00000
## GP 24.00000
## HI 28.47655
## HU 19.51429
## MS 23.31928
## NC 23.90454
## OE 21.97549
## RR 29.14014
## SC 27.53519
## WA 18.35000
## WC 24.60147structuralconnectivity<-read.csv("Reservations_LandMetrics_bbs.csv",header=T)
vegecomm_reserv_richness$reservatio## [1] "BC" "BE" "BK" "BR" "BW" "EC" "GP" "HI" "HU" "MS" "NC" "OE" "RR" "SC" "WA" "WC"structuralconnectivity$resID## [1] "BC" "BE" "BK" "BR" "BW" "EC" "GP" "HI" "HU" "MS" "NC" "OE" "RR" "SC" "WA" "WC"vegecomm_reserv_richness_struct<-cbind(vegecomm_reserv_richness,structuralconnectivity)
Rich_reserv_acres<-lm(interrich_n3~acres,vegecomm_reserv_richness_struct)
summary(Rich_reserv_acres)##
## Call:
## lm(formula = interrich_n3 ~ acres, data = vegecomm_reserv_richness_struct)
##
## Residuals:
## Min 1Q Median 3Q Max
## -3.7093 -2.6366 -0.5477 2.3302 3.5818
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 21.977011 1.013861 21.677 3.6e-12 ***
## acres 0.001382 0.000509 2.715 0.0168 *
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 2.735 on 14 degrees of freedom
## Multiple R-squared: 0.3449, Adjusted R-squared: 0.2981
## F-statistic: 7.371 on 1 and 14 DF, p-value: 0.01676Rich_reserv_ai<-lm(interrich_n3~ai,vegecomm_reserv_richness_struct)
summary(Rich_reserv_ai)##
## Call:
## lm(formula = interrich_n3 ~ ai, data = vegecomm_reserv_richness_struct)
##
## Residuals:
## Min 1Q Median 3Q Max
## -5.1392 -2.0183 0.6307 2.3067 4.5389
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 4.4278 14.9619 0.296 0.772
## ai 0.2097 0.1600 1.311 0.211
##
## Residual standard error: 3.19 on 14 degrees of freedom
## Multiple R-squared: 0.1093, Adjusted R-squared: 0.04566
## F-statistic: 1.718 on 1 and 14 DF, p-value: 0.2111Rich_reserv_areamn<-lm(interrich_n3~area_mn,vegecomm_reserv_richness_struct)
summary(Rich_reserv_areamn)##
## Call:
## lm(formula = interrich_n3 ~ area_mn, data = vegecomm_reserv_richness_struct)
##
## Residuals:
## Min 1Q Median 3Q Max
## -4.2080 -1.9539 0.2301 2.4859 4.5089
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 2.251e+01 9.114e-01 24.697 6.06e-13 ***
## area_mn 1.297e-03 5.071e-04 2.558 0.0228 *
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 2.79 on 14 degrees of freedom
## Multiple R-squared: 0.3184, Adjusted R-squared: 0.2698
## F-statistic: 6.541 on 1 and 14 DF, p-value: 0.02278Rich_reserv_ca<-lm(interrich_n3~ca,vegecomm_reserv_richness_struct)
summary(Rich_reserv_ca)##
## Call:
## lm(formula = interrich_n3 ~ ca, data = vegecomm_reserv_richness_struct)
##
## Residuals:
## Min 1Q Median 3Q Max
## -3.7135 -2.6368 -0.5517 2.3337 3.5862
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 2.198e+01 1.014e+00 21.686 3.58e-12 ***
## ca 3.169e-04 1.168e-04 2.714 0.0168 *
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 2.736 on 14 degrees of freedom
## Multiple R-squared: 0.3447, Adjusted R-squared: 0.2979
## F-statistic: 7.363 on 1 and 14 DF, p-value: 0.01681Rich_reserv_cai_mn<-lm(interrich_n3~cai_mn,vegecomm_reserv_richness_struct)
summary(Rich_reserv_cai_mn)##
## Call:
## lm(formula = interrich_n3 ~ cai_mn, data = vegecomm_reserv_richness_struct)
##
## Residuals:
## Min 1Q Median 3Q Max
## -5.8640 -1.2744 0.0397 1.4097 4.9520
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 22.01786 1.69266 13.008 3.3e-09 ***
## cai_mn 0.05042 0.03783 1.333 0.204
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 3.184 on 14 degrees of freedom
## Multiple R-squared: 0.1126, Adjusted R-squared: 0.04924
## F-statistic: 1.777 on 1 and 14 DF, p-value: 0.2038Rich_reserv_circle_mn<-lm(interrich_n3~circle_mn,vegecomm_reserv_richness_struct)
summary(Rich_reserv_circle_mn)##
## Call:
## lm(formula = interrich_n3 ~ circle_mn, data = vegecomm_reserv_richness_struct)
##
## Residuals:
## Min 1Q Median 3Q Max
## -5.2121 -1.0889 -0.1669 2.3594 4.4504
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 13.241 5.059 2.617 0.0203 *
## circle_mn 16.498 7.670 2.151 0.0494 *
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 2.93 on 14 degrees of freedom
## Multiple R-squared: 0.2484, Adjusted R-squared: 0.1947
## F-statistic: 4.627 on 1 and 14 DF, p-value: 0.04942Rich_reserv_clumpy<-lm(interrich_n3~clumpy,vegecomm_reserv_richness_struct)
summary(Rich_reserv_clumpy)##
## Call:
## lm(formula = interrich_n3 ~ clumpy, data = vegecomm_reserv_richness_struct)
##
## Residuals:
## Min 1Q Median 3Q Max
## -4.5004 -1.3690 0.5353 1.5123 3.6057
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 131.88 32.28 4.085 0.00111 **
## clumpy -105.17 31.47 -3.342 0.00484 **
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 2.521 on 14 degrees of freedom
## Multiple R-squared: 0.4438, Adjusted R-squared: 0.404
## F-statistic: 11.17 on 1 and 14 DF, p-value: 0.004841Rich_reserv_cohesion<-lm(interrich_n3~cohesion,vegecomm_reserv_richness_struct)
summary(Rich_reserv_cohesion)##
## Call:
## lm(formula = interrich_n3 ~ cohesion, data = vegecomm_reserv_richness_struct)
##
## Residuals:
## Min 1Q Median 3Q Max
## -4.7345 -1.7224 0.3427 1.9993 3.2565
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) -64.3888 30.6171 -2.103 0.0540 .
## cohesion 0.9087 0.3146 2.888 0.0119 *
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 2.675 on 14 degrees of freedom
## Multiple R-squared: 0.3733, Adjusted R-squared: 0.3286
## F-statistic: 8.34 on 1 and 14 DF, p-value: 0.01192Rich_reserv_contig_mn<-lm(interrich_n3~contig_mn,vegecomm_reserv_richness_struct)
summary(Rich_reserv_contig_mn)##
## Call:
## lm(formula = interrich_n3 ~ contig_mn, data = vegecomm_reserv_richness_struct)
##
## Residuals:
## Min 1Q Median 3Q Max
## -5.8000 -1.5313 -0.0424 1.6983 5.5598
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 20.330 2.339 8.693 5.14e-07 ***
## contig_mn 6.071 3.643 1.667 0.118
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 3.087 on 14 degrees of freedom
## Multiple R-squared: 0.1655, Adjusted R-squared: 0.1059
## F-statistic: 2.777 on 1 and 14 DF, p-value: 0.1178Rich_reserv_core_mn<-lm(interrich_n3~core_mn,vegecomm_reserv_richness_struct)
summary(Rich_reserv_core_mn)##
## Call:
## lm(formula = interrich_n3 ~ core_mn, data = vegecomm_reserv_richness_struct)
##
## Residuals:
## Min 1Q Median 3Q Max
## -4.2463 -1.9870 0.2502 2.3776 4.5108
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 2.257e+01 8.966e-01 25.17 4.67e-13 ***
## core_mn 1.470e-03 5.742e-04 2.56 0.0227 *
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 2.789 on 14 degrees of freedom
## Multiple R-squared: 0.3188, Adjusted R-squared: 0.2702
## F-statistic: 6.553 on 1 and 14 DF, p-value: 0.02268Rich_reserv_cpland<-lm(interrich_n3~cpland,vegecomm_reserv_richness_struct)
summary(Rich_reserv_cpland)##
## Call:
## lm(formula = interrich_n3 ~ cpland, data = vegecomm_reserv_richness_struct)
##
## Residuals:
## Min 1Q Median 3Q Max
## -3.8173 -2.7253 -0.4908 2.2622 3.6100
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 22.1180 0.9964 22.198 2.61e-12 ***
## cpland 0.3716 0.1411 2.634 0.0196 *
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 2.764 on 14 degrees of freedom
## Multiple R-squared: 0.3313, Adjusted R-squared: 0.2836
## F-statistic: 6.938 on 1 and 14 DF, p-value: 0.01963Rich_reserv_dcad<-lm(interrich_n3~dcad,vegecomm_reserv_richness_struct)
summary(Rich_reserv_dcad)##
## Call:
## lm(formula = interrich_n3 ~ dcad, data = vegecomm_reserv_richness_struct)
##
## Residuals:
## Min 1Q Median 3Q Max
## -4.4125 -2.8268 0.1998 1.4137 5.5563
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 22.645 1.201 18.86 2.39e-11 ***
## dcad 124.637 83.107 1.50 0.156
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 3.137 on 14 degrees of freedom
## Multiple R-squared: 0.1384, Adjusted R-squared: 0.07688
## F-statistic: 2.249 on 1 and 14 DF, p-value: 0.1559multi-species occupancy modelling
MSOM - running different models
basic null model
bbs_occobject_null.occ.formula <- ~1
bbs_occobject_null.det.formula <- ~ 1
N <- dim(bbs_occobject$y)[1]
ms.inits <- list(alpha.comm = 0,
beta.comm = 0,
beta = 0,
alpha = 0,
tau.sq.beta = 1,
tau.sq.alpha = 1,
z = apply(bbs_occobject$y, c(1, 2), max, na.rm = TRUE))
ms.priors <- list(beta.comm.normal = list(mean = 0, var = 2.72),
alpha.comm.normal = list(mean = 0, var = 2.72),
tau.sq.beta.ig = list(a = 0.1, b = 0.1),
tau.sq.alpha.ig = list(a = 0.1, b = 0.1))
out.ms.null <- msPGOcc(occ.formula = bbs_occobject_null.occ.formula,
det.formula = bbs_occobject_null.det.formula,
data = bbs_occobject,
inits = ms.inits,
n.samples = 30000,
priors = ms.priors,
n.omp.threads = 1,
verbose = TRUE,
n.report = 6000,
n.burn = 10000,
n.thin = 50,
n.chains = 3)
summary(out.ms.null)
ppc.ms.null.out <- ppcOcc(out.ms.null, 'chi-squared', group = 1)
summary(ppc.ms.null.out)
waicOcc(out.ms.null)building up complexity of detection model
bbs_occobject_null.occ.formula <- ~1
bbs_occobject_day.det.formula <- ~ scale(day)
N <- dim(bbs_occobject$y)[1]
ms.inits <- list(alpha.comm = 0,
beta.comm = 0,
beta = 0,
alpha = 0,
tau.sq.beta = 1,
tau.sq.alpha = 1,
z = apply(bbs_occobject$y, c(1, 2), max, na.rm = TRUE))
ms.priors <- list(beta.comm.normal = list(mean = 0, var = 2.72),
alpha.comm.normal = list(mean = 0, var = 2.72),
tau.sq.beta.ig = list(a = 0.1, b = 0.1),
tau.sq.alpha.ig = list(a = 0.1, b = 0.1))
out.ms.detect.day <- msPGOcc(occ.formula = bbs_occobject_null.occ.formula,
det.formula = bbs_occobject_day.det.formula,
data = bbs_occobject,
inits = ms.inits,
n.samples = 30000,
priors = ms.priors,
n.omp.threads = 1,
verbose = TRUE,
n.report = 6000,
n.burn = 10000,
n.thin = 50,
n.chains = 3)
summary(out.ms.detect.day)
ppc.ms.detect.day.out <- ppcOcc(out.ms.detect.day, 'chi-squared', group = 1)
summary(out.ms.detect.day)
waicOcc(out.ms.detect.day)day and time-of-day
bbs_occobject_null.occ.formula <- ~1
bbs_occobject_daytod.det.formula <- ~ scale(day)+scale(tod)
N <- dim(bbs_occobject$y)[1]
ms.inits <- list(alpha.comm = 0,
beta.comm = 0,
beta = 0,
alpha = 0,
tau.sq.beta = 1,
tau.sq.alpha = 1,
z = apply(bbs_occobject$y, c(1, 2), max, na.rm = TRUE))
ms.priors <- list(beta.comm.normal = list(mean = 0, var = 2.72),
alpha.comm.normal = list(mean = 0, var = 2.72),
tau.sq.beta.ig = list(a = 0.1, b = 0.1),
tau.sq.alpha.ig = list(a = 0.1, b = 0.1))
out.ms.detect.daytod <- msPGOcc(occ.formula = bbs_occobject_null.occ.formula,
det.formula = bbs_occobject_daytod.det.formula,
data = bbs_occobject,
inits = ms.inits,
n.samples = 30000,
priors = ms.priors,
n.omp.threads = 1,
verbose = TRUE,
n.report = 6000,
n.burn = 10000,
n.thin = 50,
n.chains = 3)
summary(out.ms.detect.daytod)
ppc.ms.detect.day.out <- ppcOcc(out.ms.detect.daytod, 'chi-squared', group = 1)
summary(out.ms.detect.daytod)
waicOcc(out.ms.detect.daytod)day, tod, quadratic day
bbs_occobject_null.occ.formula <- ~1
bbs_occobject_daytodquad.det.formula <- ~ scale(day)+scale(tod) +I(scale(day)^2)
N <- dim(bbs_occobject$y)[1]
ms.inits <- list(alpha.comm = 0,
beta.comm = 0,
beta = 0,
alpha = 0,
tau.sq.beta = 1,
tau.sq.alpha = 1,
z = apply(bbs_occobject$y, c(1, 2), max, na.rm = TRUE))
ms.priors <- list(beta.comm.normal = list(mean = 0, var = 2.72),
alpha.comm.normal = list(mean = 0, var = 2.72),
tau.sq.beta.ig = list(a = 0.1, b = 0.1),
tau.sq.alpha.ig = list(a = 0.1, b = 0.1))
out.ms.detect.daytodquad<- msPGOcc(occ.formula = bbs_occobject_null.occ.formula,
det.formula = bbs_occobject_daytodquad.det.formula,
data = bbs_occobject,
inits = ms.inits,
n.samples = 30000,
priors = ms.priors,
n.omp.threads = 1,
verbose = TRUE,
n.report = 6000,
n.burn = 10000,
n.thin = 50,
n.chains = 3)
summary(out.ms.detect.daytodquad)
ppc.ms.detect.day.out <- ppcOcc(out.ms.detect.daytodquad, 'chi-squared', group = 1)
summary(out.ms.detect.daytodquad)
waicOcc(out.ms.detect.daytodquad)The best-supported model is within 2 WAIC of a simpler detection model - one that just includes day. So we will proceed with just using a linear term for day in our detectability submodels.
Now we will run several different occupancy submodels, once again evaluating fit for progressively more complicated models.
occupancy model with vibi-fq
bbs_occobject_vibifq.occ.formula <- ~ scale(vibifq)
bbs_occobject.det.best.formula <- ~ scale(day)
N <- dim(bbs_occobject$y)[1]
ms.inits <- list(alpha.comm = 0,
beta.comm = 0,
beta = 0,
alpha = 0,
tau.sq.beta = 1,
tau.sq.alpha = 1,
z = apply(bbs_occobject$y, c(1, 2), max, na.rm = TRUE))
ms.priors <- list(beta.comm.normal = list(mean = 0, var = 2.72),
alpha.comm.normal = list(mean = 0, var = 2.72),
tau.sq.beta.ig = list(a = 0.1, b = 0.1),
tau.sq.alpha.ig = list(a = 0.1, b = 0.1))
out.ms.vibifq <- msPGOcc(occ.formula = bbs_occobject_vibifq.occ.formula,
det.formula = bbs_occobject.det.best.formula,
data = bbs_occobject,
inits = ms.inits,
n.samples = 30000,
priors = ms.priors,
n.omp.threads = 1,
verbose = TRUE,
n.report = 6000,
n.burn = 10000,
n.thin = 50,
n.chains = 3)
summary(out.ms.vibifq)
ppc.ms.vibifq.out <- ppcOcc(out.ms.vibifq, 'chi-squared', group = 1)
summary(ppc.ms.vibifq.out)
out.ms.bestdec.small <- msPGOcc(occ.formula = ~1,
det.formula = bbs_occobject.det.best.formula,
data = bbs_occobject,
inits = ms.inits,
n.samples = 30000,
priors = ms.priors,
n.omp.threads = 1,
verbose = TRUE,
n.report = 6000,
n.burn = 10000,
n.thin = 50,
n.chains = 3)
waicOcc(out.ms.bestdec.small)
waicOcc(out.ms.vibifq)occupancy model with vibi-fq-notrees
bbs_occobject_vibifqnotrees.occ.formula <- ~ scale(vibifqnotrees)
bbs_occobject.det.best.formula <- ~ scale(day)
N <- dim(bbs_occobject$y)[1]
ms.inits <- list(alpha.comm = 0,
beta.comm = 0,
beta = 0,
alpha = 0,
tau.sq.beta = 1,
tau.sq.alpha = 1,
z = apply(bbs_occobject$y, c(1, 2), max, na.rm = TRUE))
ms.priors <- list(beta.comm.normal = list(mean = 0, var = 2.72),
alpha.comm.normal = list(mean = 0, var = 2.72),
tau.sq.beta.ig = list(a = 0.1, b = 0.1),
tau.sq.alpha.ig = list(a = 0.1, b = 0.1))
out.ms.vibifqnotrees <- msPGOcc(occ.formula = bbs_occobject_vibifqnotrees.occ.formula,
det.formula = bbs_occobject.det.best.formula,
data = bbs_occobject,
inits = ms.inits,
n.samples = 30000,
priors = ms.priors,
n.omp.threads = 1,
verbose = TRUE,
n.report = 6000,
n.burn = 10000,
n.thin = 50,
n.chains = 3)
summary(out.ms.vibifqnotrees)
ppc.ms.vibifqnotrees.out <- ppcOcc(out.ms.vibifqnotrees, 'chi-squared', group = 1)
summary(ppc.ms.vibifqnotrees.out)
waicOcc(out.ms.vibifqnotrees)occupancy model with all rda-supported covariates
bbs_occobject_fullrda.occ.formula <- ~ scale(tol_rel) + scale(svp) + scale(dicot) + scale(cyper) + scale(small_tree)
bbs_occobject.det.best.formula <- ~ scale(day)
N <- dim(bbs_occobject$y)[1]
ms.inits <- list(alpha.comm = 0,
beta.comm = 0,
beta = 0,
alpha = 0,
tau.sq.beta = 1,
tau.sq.alpha = 1,
z = apply(bbs_occobject$y, c(1, 2), max, na.rm = TRUE))
ms.priors <- list(beta.comm.normal = list(mean = 0, var = 2.72),
alpha.comm.normal = list(mean = 0, var = 2.72),
tau.sq.beta.ig = list(a = 0.1, b = 0.1),
tau.sq.alpha.ig = list(a = 0.1, b = 0.1))
out.ms.fullrda <- msPGOcc(occ.formula = bbs_occobject_fullrda.occ.formula,
det.formula = bbs_occobject.det.best.formula,
data = bbs_occobject,
inits = ms.inits,
n.samples = 30000,
priors = ms.priors,
n.omp.threads = 1,
verbose = TRUE,
n.report = 6000,
n.burn = 10000,
n.thin = 50,
n.chains = 3)
summary(out.ms.fullrda)
ppc.ms.fullrda.out <- ppcOcc(out.ms.fullrda, 'chi-squared', group = 1)
summary(ppc.ms.fullrda.out)
waicOcc(out.ms.fullrda)occupancy model with just tolerance/sensitivity
bbs_occobject_tolerance.occ.formula <- ~ scale(tol_rel)
bbs_occobject.det.best.formula <- ~ scale(day)
N <- dim(bbs_occobject$y)[1]
ms.inits <- list(alpha.comm = 0,
beta.comm = 0,
beta = 0,
alpha = 0,
tau.sq.beta = 1,
tau.sq.alpha = 1,
z = apply(bbs_occobject$y, c(1, 2), max, na.rm = TRUE))
ms.priors <- list(beta.comm.normal = list(mean = 0, var = 2.72),
alpha.comm.normal = list(mean = 0, var = 2.72),
tau.sq.beta.ig = list(a = 0.1, b = 0.1),
tau.sq.alpha.ig = list(a = 0.1, b = 0.1))
out.ms.tolerance<- msPGOcc(occ.formula = bbs_occobject_tolerance.occ.formula,
det.formula = bbs_occobject.det.best.formula,
data = bbs_occobject,
inits = ms.inits,
n.samples = 30000,
priors = ms.priors,
n.omp.threads = 1,
verbose = TRUE,
n.report = 6000,
n.burn = 10000,
n.thin = 50,
n.chains = 3)
summary(out.ms.tolerance)
ppc.ms.tolerance.out <- ppcOcc(out.ms.tolerance, 'chi-squared', group = 1)
summary(ppc.ms.tolerance.out)
waicOcc(out.ms.tolerance)global model
bbs_occobject_global.occ.formula <- ~ scale(tol_rel) + scale(svp) + scale(dicot) + scale(cyper) + scale(small_tree)+scale(shrub)+scale(hydro)+scale(subcanopy)+scale(canopy)
bbs_occobject.det.best.formula <- ~ scale(day)
N <- dim(bbs_occobject$y)[1]
ms.inits <- list(alpha.comm = 0,
beta.comm = 0,
beta = 0,
alpha = 0,
tau.sq.beta = 1,
tau.sq.alpha = 1,
z = apply(bbs_occobject$y, c(1, 2), max, na.rm = TRUE))
ms.priors <- list(beta.comm.normal = list(mean = 0, var = 2.72),
alpha.comm.normal = list(mean = 0, var = 2.72),
tau.sq.beta.ig = list(a = 0.1, b = 0.1),
tau.sq.alpha.ig = list(a = 0.1, b = 0.1))
out.ms.global <- msPGOcc(occ.formula = bbs_occobject_global.occ.formula,
det.formula = bbs_occobject.det.best.formula,
data = bbs_occobject,
inits = ms.inits,
n.samples = 30000,
priors = ms.priors,
n.omp.threads = 1,
verbose = TRUE,
n.report = 6000,
n.burn = 10000,
n.thin = 50,
n.chains = 3)
summary(out.ms.global)
ppc.ms.global.out <- ppcOcc(out.ms.global, 'chi-squared', group = 1)
summary(ppc.ms.global.out)
waicOcc(out.ms.global)just FQAI
bbs_occobject_fqi.occ.formula <- ~ scale(fqi)
bbs_occobject.det.best.formula <- ~ scale(day)
N <- dim(bbs_occobject$y)[1]
ms.inits <- list(alpha.comm = 0,
beta.comm = 0,
beta = 0,
alpha = 0,
tau.sq.beta = 1,
tau.sq.alpha = 1,
z = apply(bbs_occobject$y, c(1, 2), max, na.rm = TRUE))
ms.priors <- list(beta.comm.normal = list(mean = 0, var = 2.72),
alpha.comm.normal = list(mean = 0, var = 2.72),
tau.sq.beta.ig = list(a = 0.1, b = 0.1),
tau.sq.alpha.ig = list(a = 0.1, b = 0.1))
out.ms.fqi <- msPGOcc(occ.formula = bbs_occobject_fqi.occ.formula,
det.formula = bbs_occobject.det.best.formula,
data = bbs_occobject,
inits = ms.inits,
n.samples = 30000,
priors = ms.priors,
n.omp.threads = 1,
verbose = TRUE,
n.report = 6000,
n.burn = 10000,
n.thin = 50,
n.chains = 3)
summary(out.ms.fqi)
ppc.ms.fqi.out <- ppcOcc(out.ms.fqi, 'chi-squared', group = 1)
summary(ppc.ms.fqi.out)
waicOcc(out.ms.fqi)model selection
we can also print out a table of Widely Applicable Information Criterion (wAIC) values. The lower this is, the better the model fit.
AICc_table<-tibble(.rows = 10)
AICc_table$model<-c("Null","p ~ day", "p ~ day + time-of-day", "p ~ day^2 + day +tod", "p ~ day, psi ~ vibi-fq", "p ~ day, psi ~ vibi-fq-notree","p ~ day, psi ~ tolerant + seedless vascular + dicots + cyperaceae + small trees","p~day,psi~tolerant","p~day,psi~tolerant + seedless vascular + dicot + cyperaceae + small trees +shrubs+hydro+subcanopy+canopy","p~day,psi~fqi")
AICc_table$wAIC<-c(waicOcc(out.ms.null)[[3]],waicOcc(out.ms.detect.day)[[3]],waicOcc(out.ms.detect.daytod)[[3]],
waicOcc(out.ms.detect.daytodquad)[[3]],waicOcc(out.ms.vibifq)[[3]],waicOcc(out.ms.vibifqnotrees)[[3]],
waicOcc(out.ms.fullrda)[[3]],waicOcc(out.ms.tolerance)[[3]],waicOcc(out.ms.global)[[3]],waicOcc(out.ms.fqi))## X model wAIC
## 9 9 p~day,psi~tolerant + seedless vascular + dicot + cyperaceae + small trees +shrubs+hydro+subcanopy+canopy 23165.73
## 7 7 p ~ day, psi ~ tolerant + seedless vascular + dicots + cyperaceae + small trees 23185.00
## 5 5 p ~ day, psi ~ vibi-fq 23190.57
## 8 8 p~day,psi~tolerant 23323.25
## 6 6 p ~ day, psi ~ vibi-fq-notree 23372.44
## 4 4 p ~ day^2 + day +tod 23755.71
## 2 2 p ~ day 23756.93
## 3 3 p ~ day + time-of-day 23769.78
## 1 1 Null 23840.32So, the top model is the global, followed by the fullrda model.
Now, we plot based on this top model.
plotting using an RShiny app
covariate plotting
#options(max.print=1000000)
#sink("outmsglobal_90percent.txt")
#print(summary(out.ms.global,level="species",quantiles=c(0.05,0.5,0.95)))
#sink()
covars<-read.csv("outmsglobal_covar.csv",header=T)
covars_tolerance<-subset(covars,covars$Covariate=="tol_rel")
covars_tolerance<-covars_tolerance[order(covars_tolerance$Mean),]
spp_tolorder<-factor(covars_tolerance$Species,levels=covars_tolerance$Species)
covars_tolerance$spporder<-spp_tolorder
ggplot(aes(x=Mean,y=spp_tolorder),data=covars_tolerance)+geom_point()+geom_errorbarh(aes(xmin=X5.0CI,xmax=X95.0CI))+geom_vline(xintercept=0.1652)+theme_classic(base_size=16)
covars_svp<-subset(covars,covars$Covariate=="svp")
covars_svp<-covars_svp[order(covars_svp$Mean),]
spp_svporder<-factor(covars_svp$Species,levels=covars_svp$Species)
covars_svp$spporder<-spp_svporder
ggplot(aes(x=Mean,y=spp_svporder),data=covars_svp)+geom_point()+geom_errorbarh(aes(xmin=X5.0CI,xmax=X95.0CI))+geom_vline(xintercept=-0.0594)+theme_classic(base_size=16)
covars_dicot<-subset(covars,covars$Covariate=="dicot")
covars_dicot<-covars_dicot[order(covars_dicot$Mean),]
spp_dicotorder<-factor(covars_dicot$Species,levels=covars_dicot$Species)
covars_dicot$spporder<-spp_dicotorder
ggplot(aes(x=Mean,y=spp_dicotorder),data=covars_dicot)+geom_point()+geom_errorbarh(aes(xmin=X5.0CI,xmax=X95.0CI))+geom_vline(xintercept=0.2461)+theme_classic(base_size=16)
covars_cyper<-subset(covars,covars$Covariate=="cyper")
covars_cyper<-covars_cyper[order(covars_cyper$Mean),]
spp_cyperorder<-factor(covars_cyper$Species,levels=covars_cyper$Species)
covars_cyper$spporder<-spp_cyperorder
ggplot(aes(x=Mean,y=spp_cyperorder),data=covars_cyper)+geom_point()+geom_errorbarh(aes(xmin=X5.0CI,xmax=X95.0CI))+geom_vline(xintercept=-0.1375)+theme_classic(base_size=16)
covars_smalltree<-subset(covars,covars$Covariate=="small_tree")
covars_smalltree<-covars_smalltree[order(covars_smalltree$Mean),]
spp_treeorder<-factor(covars_smalltree$Species,levels=covars_smalltree$Species)
covars_smalltree$spporder<-spp_treeorder
ggplot(aes(x=Mean,y=spp_treeorder),data=covars_smalltree)+geom_point()+geom_errorbarh(aes(xmin=X5.0CI,xmax=X95.0CI))+geom_vline(xintercept=0.0672)+theme_classic(base_size=16)
covars_shrub<-subset(covars,covars$Covariate=="shrub")
covars_shrub<-covars_shrub[order(covars_shrub$Mean),]
spp_shruborder<-factor(covars_shrub$Species,levels=covars_shrub$Species)
covars_shrub$spporder<-spp_shruborder
ggplot(aes(x=Mean,y=spp_shruborder),data=covars_shrub)+geom_point()+geom_errorbarh(aes(xmin=X5.0CI,xmax=X95.0CI))+geom_vline(xintercept=-0.2181)+theme_classic(base_size=16)
covars_hydro<-subset(covars,covars$Covariate=="hydro")
covars_hydro<-covars_hydro[order(covars_hydro$Mean),]
spp_hydroorder<-factor(covars_hydro$Species,levels=covars_hydro$Species)
covars_hydro$spporder<-spp_hydroorder
ggplot(aes(x=Mean,y=spp_hydroorder),data=covars_hydro)+geom_point()+geom_errorbarh(aes(xmin=X5.0CI,xmax=X95.0CI))+geom_vline(xintercept=0.2068)+theme_classic(base_size=16)
covars_subcanopy<-subset(covars,covars$Covariate=="subcanopy")
covars_subcanopy<-covars_subcanopy[order(covars_subcanopy$Mean),]
spp_subcanopyorder<-factor(covars_subcanopy$Species,levels=covars_subcanopy$Species)
covars_subcanopy$spporder<-spp_subcanopyorder
ggplot(aes(x=Mean,y=spp_subcanopyorder),data=covars_subcanopy)+geom_point()+geom_errorbarh(aes(xmin=X5.0CI,xmax=X95.0CI))+geom_vline(xintercept=-0.0393)+theme_classic(base_size=16)
covars_canopy<-subset(covars,covars$Covariate=="canopy")
covars_canopy<-covars_canopy[order(covars_canopy$Mean),]
spp_canopyorder<-factor(covars_canopy$Species,levels=covars_canopy$Species)
covars_canopy$spporder<-spp_canopyorder
ggplot(aes(x=Mean,y=spp_canopyorder),data=covars_canopy)+geom_point()+geom_errorbarh(aes(xmin=X5.0CI,xmax=X95.0CI))+geom_vline(xintercept=-0.0702)+theme_classic(base_size=16)